Filters look good - have you tried verifying by convoluting speakerrespone with correction filter and comparing it with you desired response (all in time domain)??
I'll have a look on a low frequency filter soon - or if you want to do it yourself; here is the code:
clf(1)
//Definitions
ntaps = 500
fs = 44100
freqvect = [0:ntaps-1]*fs/2/ntaps
//Import measured data and the desired response (target)
//Measure data is assumed to be three columns of data: frequency (Hz), magnitude (dB), phase (rad).
//For the target however, I choose to have a linear scale for the magnitude, since I'd like zero for the stopband (which gets -inf in dB)
measureddata = fscanfMat('/path/measured.txt');
targetdata = fscanfMat('/path/target.txt');
//For calculation of the "ideal" filter - the imported target data is interpolated to the same frequency vector as the measured data
f=measureddata(:,1) //The frequency vector
//Reform the measured data into complex rectangular form
magn = measureddata(:,2)
phase = measureddata(:,3)
Repart = 10^(magn/20).*cos(phase*2*%pi/360)
Impart = 10^(magn/20).*sin(phase*2*%pi/360)
measured_complex = [Repart+Impart*%i];
//Interpolate the target function for the same frequency
df = splin(targetdata(:,1),targetdata(:,2))
dp = splin(targetdata(:,1),targetdata(:,3))
magn = interp(f,targetdata(:,1),targetdata(:,2),df)
phase = interp(f, targetdata(:,1),targetdata(:,3),dp)
Repart = magn.*cos(phase*2*%pi/360)
Impart = magn.*sin(phase*2*%pi/360)
target_complex = [Repart+Impart*%i];
//Now, create the "ideal" filter by dividing the target with the measured data
ideal_complex = target_complex./measured_complex
//The ideal filter is truncated in the time domain by resampling in the frequency domain. The interp function require real numbers, so the complex ideal filter must first be separated into real and imaginary parts
Re=real(ideal_complex)
Im=imag(ideal_complex)
dr = splin(f,Re)
di = splin(f,Im)
Re_trunc = interp(freqvect,f,Re,dr)
Im_trunc = interp(freqvect,f,Im,di)
//Combine the resampled real vectors into the new complex filter response vector
filter_trunc = Re_trunc+Im_trunc*%i
//Do ifft on the filter response and shift the result
filter_trunc_ifft = (fftshift(ifft(filter_trunc)))
//Apply the window to the truncated filter (in this case a Hamming window)
wind_filter=filter_trunc_ifft.*window('kr',ntaps,5)
//Transform back to get the frequency and phase response of the filter
filter_response = fft(wind_filter)
//Calculating the magnitude response of the filter
filter_magn = abs(filter_response(1:ntaps))
//Calculating the magnitude response of the "ideal" filter
ideal_response = abs(ideal_complex)
// If there are any zero elements the logarithm can't be calculated, so zeros are replaced by a very small number
x=find (ideal_response<=0.) ; ideal_response(x)=0.00000001
x=find (targetdata(:,2)<=0.) ; targetdata(x,2)=0.00000001
//Calculate the resulting response after filtering by convolution in the time domain
Result = fft(convol((wind_filter),fftshift(ifft(measured_complex))))
//Plot the result
scf(1)
subplot(2,3,1)
plot2d ('ln',measureddata(:,1),measureddata(:,2),rect=[20.,-100.,20000.,10.]);title("Imported measured response",'fontsize',3)
subplot(2,3,2)
plot2d ('ln',targetdata(:,1),20*log10(targetdata(:,2)),rect=[20.,-100.,20000.,10.]);title("Imported target response",'fontsize',3)
subplot(2,3,3)
plot2d ('ln',f,20*log10(ideal_response),rect=[20.,-100.,20000.,10.]);title("Ideal filter response",'fontsize',3)
subplot(2,3,4)
plot2d(wind_filter);title("Truncated filter",'fontsize',3)
subplot(2,3,5)
plot2d ('ln',freqvect,20*log10(filter_magn),rect=[20.,-100.,20000.,10.]);title("Resulting filter response",'fontsize',3)
subplot(2,3,6)
plot2d ('ln',f,20*log10(abs(Result(1:length(f)))),rect=[20.,-100.,20000.,10.]);title("Result of convolution measured data with filter",'fontsize',3)