progress on IIR

matlab/ITADSPThirdOctaveFilterbankCoefficients.proto.h

+/*
+ * ----------------------------------------------------------------
+ *
+ *		ITA core libs
+ *		(c) Copyright Institute of Technical Acoustics (ITA)
+ *		RWTH Aachen University, Germany, 2015-2017
+ *
+ * ----------------------------------------------------------------
+ *				    ____  __________  _______
+ *				   //  / //__   ___/ //  _   |
+ *				  //  /    //  /    //  /_|  |
+ *				 //  /    //  /    //  ___   |
+ *				//__/    //__/    //__/   |__|
+ *
+ * ----------------------------------------------------------------
+ *
+ */
+
+### stub code ###

matlab/ita_dsp_design_filters.m

-
+%% Settings
 Fs = 44100;
-
-Hd = ita_dsp_third_octave_filterbank_iir_lowpass_coefficients( Fs, 14.1 )
-ita_dsp_third_octave_filterbank_iir_band_coefficients( Fs, 14.1, 17.8 ) % 16 Hz band
-% ...
-ita_dsp_third_octave_filterbank_iir_highpass_coefficients( Fs, 17780 )
+N = 10;
+freqs = ita_ANSI_center_frequencies;
 
 %% Export to C/C++ header using std::vector
 
-copyfile( 'ITADSPThirdOctaveFilterbankCoefficients.proto.h', 'ITADSPThirdOctaveFilterbankCoefficients.h' );
-fd = fopen( 'ITADSPThirdOctaveFilterbankCoefficients.h', 'wt' );
-fseek( fd, '### stub code ###' );
-fprintf( fd, 'std::vector< float > vfITADSPThirdOctaveFilterbankCoefficients = {' );
-fprintf( fs, '%.9g', length( Hd.sosMatrix ), h( 1 ) );
-fprintf( fd, ',\n %.9g', h( 2:end ) );
-fprintf( fd, '};\n');
+target_file = 'ITADSPThirdOctaveFilterbankIIRCoefficients.h';
+proto_file = 'ITADSPThirdOctaveFilterbankCoefficients.proto.h';
+copyfile( proto_file, target_file );
+
+fd = fopen( target_file, 'r' );
+stub_cpp_file_body = textscan( fd, '%s', 'delimiter', '\n' );
+fclose( fd );
+
+fd = fopen( target_file, 'wt' );
+for i=1:size( stub_cpp_file_body{ 1 } , 1 )
+    if strcmp( '### stub code ###',  stub_cpp_file_body{ 1 }{ i } ) 
+        break;
+    end
+    fprintf( fd, '%s\n', stub_cpp_file_body{ 1 }{ i } );
+end
+
+fprintf( fd, '/* This file has been generated by the ita_dsp_design_filsters.m Matlab script */\n\n' );
+
+fprintf( fd, '#ifndef IW_ITA_DSP_THIRD_OCTAVE_FILTERBANK_IIR_COEFFICIENTS\n' );
+fprintf( fd, '#define IW_ITA_DSP_THIRD_OCTAVE_FILTERBANK_IIR_COEFFICIENTS\n\n' );
+fprintf( fd, '#include <vector>\n\n' );
+
+fprintf( fd, 'namespace ITADSPThirdOctaveFilterbankIIR\n{\n\n' );
+fprintf( fd, 'const double dSamplingFrequency = %0.2f; // Hz\n', Fs );
+fprintf( fd, 'const int iOrdner = %d;\n', N );
+fprintf( fd, 'const int iNumFrequencies = %d;\n', numel( freqs ) );
+fprintf( fd, 'const std::vector< float > vfFrequencies =\n{\n' );
+for l=1:numel( freqs )
+    fprintf( fd, '%.1ff, // Hz\n', freqs( l ) );
+end
+fprintf( fd, '};\n' );
+fprintf( fd, '\n' );
+
+%% Filter design
+
+fprintf( fd, '// 3-dim vector access: FrequencyBand ; Biquad ; Param = { ScaleValue [1], sosMatrix [2-6] }\n' );
+fprintf( fd, 'const std::vector< std::vector< std::vector< float > > > vvvfSOSMatrixCoefficients = \n{\n\n' );
+
+% Lowpass
+h  = fdesign.lowpass( 'N,F3dB', N, freqs( 1 ), Fs );
+Hd = design( h, 'butter' );
+
+fprintf( fd, '{\n' );
+fprintf( fd, '// Band %.1f Hz sos matrix\n', freqs( 1 ) );
+for k=1:size( Hd.sosMatrix, 1 )
+    fprintf( fd, '{ %.6ff, ', Hd.ScaleValues( k ) );
+    for l=1:size( Hd.sosMatrix, 2 )
+        fprintf( fd, '%.6ff, ', Hd.sosMatrix( k, l ) );
+    end
+    fprintf( fd, '},\n' );
+end
+fprintf( fd, '},\n\n' );
+
+% Bandpasses
+for n = 2:( numel( freqs ) - 1 )
+    
+    h  = fdesign.bandpass( 'N,F3dB1,F3dB2', N, freqs( n - 1 ), freqs( n + 1 ), Fs );
+    Hd = design( h, 'butter' );
+
+    fprintf( fd, '{\n' );
+    fprintf( fd, '// Band %.1f Hz sos matrix\n', freqs( n ) );
+    for k=1:size( Hd.sosMatrix, 1 )
+        fprintf( fd, '{ %.6ff, ', Hd.ScaleValues( k ) );
+        for l=1:size( Hd.sosMatrix, 2 )
+            fprintf( fd, '%.6ff, ', Hd.sosMatrix( k, l ) );
+        end
+        fprintf( fd, '},\n' );
+    end
+    fprintf( fd, '},\n\n' );
+    
+end
+
+% Highpass
+h  = fdesign.highpass( 'N,F3dB', N, freqs( end ), Fs );
+Hd = design( h, 'butter' );
+
+fprintf( fd, '{\n' );
+fprintf( fd, '// Band %.1f Hz sos matrix\n', freqs( end ) );
+for k=1:size( Hd.sosMatrix, 1 )
+	fprintf( fd, '{ %.6ff, ', Hd.ScaleValues( k ) );
+    for l=1:size( Hd.sosMatrix, 2 )
+        fprintf( fd, '%.6ff, ', Hd.sosMatrix( k, l ) );
+    end
+    fprintf( fd, '},\n' );
+end
+fprintf( fd, '},\n\n' );
+
+fprintf( fd, '};\n\n' );
+
+fprintf( fd, '} // namespace\n\n' );
+
+fprintf( fd, '#endif // IW_ITA_DSP_THIRD_OCTAVE_FILTERBANK_IIR_COEFFICIENTS\n' );
+
+%% Finish
+for ii=(i+1):size( stub_cpp_file_body{ 1 } , 1 )
+    fprintf( fd, '%s\n', stub_cpp_file_body{ 1 }{ ii } );
+end
 fclose( fd );

matlab/ita_dsp_third_octave_filterbank_fir_band_coefficients.m

+function Hd = ita_dsp_third_octave_filterbank_fir_band_coefficients( Fs, Fc1, Fc2 )
+%ITA_DSP_THIRD_OCTAVE_FILTERBANK_IIR_COEFFICIENTS Returns a discrete-time filter object.
+%   Fs Sampling Frequency
+%   Fc1 First Cutoff Frequency
+%   Fc2 Second Cutoff Frequency
+
+N   = 10;     % Order
+d = fdesign.bandpass( 'Fst1,Fp1,Fp2,Fst2,Ast1,Ap,Ast2', 1/4, 3/8, 5/8, 6/8, 60, 1, 60, Fs );
+Hd = design(d,'equiripple');