diff options
| author | David Doukhan <rudk@users.sourceforge.net> | 2011-02-07 13:09:27 +0000 |
|---|---|---|
| committer | David Doukhan <rudk@users.sourceforge.net> | 2011-02-07 13:09:27 +0000 |
| commit | 45dfb5debd4a0d0bc79ad8b6691576e8462d1396 (patch) | |
| tree | 7b4d1736586da6ac25a87d5135a357bab291648e | |
| parent | d69fd2f6b4cbbc538113e107c5d52a5d27d197bd (diff) | |
cw_binaural~ code and examples!
Makefile still missing
svn path=/trunk/externals/ddoukhan/cw_binaural~/; revision=14855
44 files changed, 3057 insertions, 0 deletions
diff --git a/README.txt b/README.txt new file mode 100644 index 0000000..611feeb --- /dev/null +++ b/README.txt @@ -0,0 +1,14 @@ +cw_binaural~: A binaural synthesis external + +Authors: + David Doukhan <david.doukhan@gmail.com> - http://www.limsi.fr/Individu/doukhan + Anne Sedes <sedes.anne@gmail.com> + +For more details, see: +CW_binaural: a binaural Sythesis External for Pure Data +David Doukhan and Anne Sedes, PDCON09 + +* libsndfile should be installed on your machine in order to compile the external + +* The external use try/catch statements. On some OS, it may cause PD to crash when exceptions are thrown. To solve this issue, go to file=>path=>startup, edit a new entry, and set cw_binaural~ to be the first external to be loaded. +For more details concerning the try/cath issue, see http://www.mail-archive.com/pd-dev@iem.at/msg06694.html diff --git a/examples/01-getting_started.pd b/examples/01-getting_started.pd new file mode 100644 index 0000000..6154fd3 --- /dev/null +++ b/examples/01-getting_started.pd @@ -0,0 +1,35 @@ +#N canvas 277 247 461 596 10; +#X obj 122 299 cw_binaural~; +#X obj 144 353 dac~; +#X obj 141 245 *~ 360; +#X obj 37 235 noise~; +#X obj 133 221 phasor~ 0.1; +#X floatatom 290 227 5 0 0 0 - - -; +#X msg 32 430 listen_db /home/david/listen/1048/COMPENSATED/WAV/IRC_1048_C +; +#X text 19 212 audio input; +#X text 132 201 azimuth (degree); +#X text 271 201 elevation (degree); +#X text 67 387 A click on this message loads the HRTF that will be +used: it is required to perform any processing.; +#X text 44 459 Listen HRIRs can be downloaded on the website of IRCAM's +Room Acoustics Team: see http://recherche.ircam.fr/equipes/salles/listen +Once a given HRIR has been unzipped \, the path to the directory containing +wavs corresponding to conpensated hrir should be given in the message. +; +#X obj 178 119 image azimuth.gif; +#X obj 350 133 image elevation.gif; +#X text 53 27 azimuth and elevation are expressed in vertical polar +coordinates; +#X text 59 -61 This patch performs binaural synthesis using cw_binaural~ +; +#X text 58 -37 More details concerning the methods used for the synthesis +can be found in: cw_binaural~: A binaural synthesis external for pure +data \, David Doukhan and Anne Sedes \, PDCon09; +#X connect 0 0 1 0; +#X connect 0 1 1 1; +#X connect 2 0 0 1; +#X connect 3 0 0 0; +#X connect 4 0 2 0; +#X connect 5 0 0 2; +#X connect 6 0 0 0; diff --git a/examples/02-supported_hrtf_database.pd b/examples/02-supported_hrtf_database.pd new file mode 100644 index 0000000..412375e --- /dev/null +++ b/examples/02-supported_hrtf_database.pd @@ -0,0 +1,23 @@ +#N canvas 1108 244 438 547 10; +#X text 60 9 Support for HRTF Databases; +#X msg -18 278 listen_db /yourpath/1048/COMPENSATED/WAV/IRC_1048_C +; +#X text -15 125 Listen HRTF Database it can be downloaded on the website +of IRCAM's Room Acoustics Team: see http://recherche.ircam.fr/equipes/salles/listen +Once a given HRIR has been unzipped \, the path to the directory containing +wavs corresponding to conpensated hrir should be given.; +#X text -21 317 CIPIC HRTF Database the database is available in .mat +format at http://interface.idav.ucdavis.edu/data/CIPIC_hrtf_database.zip +it can be converted to wav files using the matlab script cipicdb2wav.m +provided in cw_binaural archive Otherwise \, you can download the converted +database at http://www.limsi.fr/Individu/doukhan/cipic.bz2; +#X msg -21 490 cipic_db /yourpath/subject_003; +#X text -3 43 cw_binaural~ provide suppport for any hrtf database \, +assuming each impulse response is stored in a wav file whose name contains +the azimuth and elevation of the impulse response in degree.; +#X text -14 223 the syntax of the message to send to the external to +load a listen hrtf database is: 'listen_db' followed by the path to +the compensated impulse response; +#X text -23 439 the syntax of the message to send to the external to +load a CIPIC hrtf database is: 'cipic_db' followed by the path to the +compensated impulse response; diff --git a/examples/03-advanded_support_of_hrtf_database.pd b/examples/03-advanded_support_of_hrtf_database.pd new file mode 100644 index 0000000..94d83aa --- /dev/null +++ b/examples/03-advanded_support_of_hrtf_database.pd @@ -0,0 +1,29 @@ +#N canvas 175 52 512 755 10; +#X text 100 16 Advanded Support for HRTF database; +#X obj 222 408 image coordinates.gif; +#X msg 19 549 listen_db /yourpath/1048/COMPENSATED/WAV/IRC_1048_C; +#X text 14 527 the two messages bellow are equivalent to load listen +HRTF; +#X text 24 637 the two messages bellow are equivalent to load CIPIC +HRTF; +#X msg 20 572 set_hrtf_db /yourpath/1048/COMPENSATED/WAV/IRC_1048_C +IRC_[0-9]+_C_R[0-9]+_T([-+]?[0-9]*//.?[0-9]+)_P([-+]?[0-9]*//.?[0-9]+)//.wav +true true; +#X msg 22 661 cipic_db /yourpath/subject_003; +#X msg 22 688 set_hrtf_db /yourpath/subject_003 subject_[0-9]+_azim_([-+]?[0-9]//*.?[0-9]+)_elev_([-+]?[0-9]*//.?[0-9]+)//.wav +true false; +#X text 20 46 cw_binaural~ provide support for any hrtf database that +can be used through the message set_hrtf_db \, which has 4 arguments +1st arg is the path to the directory containing the hrtf database to +use \, such than each impulse response is stored in a wave file whose +name contains the azimuth and elevation of the hrtf in degree 2nd arg +is a regexp telling how to decode the wav file name: the regexp MUST +have 2 groups corresponding to the azimuth and elevation 3rd arg should +have the values 'true' \, or 'false': it tells if the azimuth correspond +to the first group of the regexp. if set to false \, then the elevation +is the first group of the regexp 4th arg should have the values 'true' +or 'false'. If set to true \, then the database is sampled using vertical-polar +coordinates (like listen database) \, if set to false \, the database +is sampled using interaural-polar coordinates (like cipic database). +Those coordinate systems are represented in the figure bellow \, comming +from the CIPIC Interface lab.; diff --git a/examples/04-using_subset_of_hrtf_database.pd b/examples/04-using_subset_of_hrtf_database.pd new file mode 100644 index 0000000..1bb978c --- /dev/null +++ b/examples/04-using_subset_of_hrtf_database.pd @@ -0,0 +1,8 @@ +#N canvas 436 209 421 196 10; +#X text 103 14 Using subset of HRTF database; +#X text 42 60 When used for research aims \, it may be usefull to consider +only a subset of the available impulse responses. In that aim \, you +just need to make a directory containing the wav files corersponding +to the subset of the impulse responses you want to use \, and give +the path of that directory through the 'listen_db' \, 'cipic_dc' or +'set_hrtf_db_messages'; diff --git a/examples/05-cw_binaural~_creation_arguments.pd b/examples/05-cw_binaural~_creation_arguments.pd new file mode 100644 index 0000000..2ab4bba --- /dev/null +++ b/examples/05-cw_binaural~_creation_arguments.pd @@ -0,0 +1,24 @@ +#N canvas 855 359 391 423 10; +#X text 105 18 cw_binaural~ creation arguments; +#X text 21 50 each instance of the external has 3 creation arguments +with defaults values; +#X text 21 83 1st arg: is the length of the impulse response used for +filtering \, default value is 128; +#X text 22 115 2nd arg: is the filtering method: default value is 'RIFF' +\, which means that the filtering will be done in temporal domain: +this method has a huge computational cost \, but 0 ltency. The other +supported filtering method is 'FFT' \, it lowers the filtering cost +\, implies more latency \, and require to use the external in a block +\, examples will be given later.; +#X text 23 198 3rd arg: if set to 'nodelay' then interpolation is done +directly on the HRIR \, and results in less computations \, but results +in a not so good itd estimation in the interpolated filter. Otherwise +\, we proceed to a minimalphase and pure delay decomposition of the +HRIR. The user may chose the fractional delay algorithm used between +'Hermite4' \, '6points' \, 'linear' \, and 'nofractional'. Default +value for the 3rd arg is 'Hermite4' \, '6points' is also a good method +but more expensive. 'linear' and 'nofractional' should not be used +\, except for didactic purposes.; +#X obj 74 359 cw_binaural~ 128 RIFF Hermite4; +#X text 36 378 instantiation of the external with its defaults creation +arguments; diff --git a/examples/06-FFT_filtering.pd b/examples/06-FFT_filtering.pd new file mode 100644 index 0000000..0bbfde3 --- /dev/null +++ b/examples/06-FFT_filtering.pd @@ -0,0 +1,20 @@ +#N canvas 350 242 427 342 10; +#X text 73 13 doing the filtering in the spectral domain; +#X text 36 49 the filtering is the costliest step of the computations +and can be done in the spectral domain to lower the amound of ressources +neccessary. The patch below show how to use it.; +#X obj 134 238 dac~; +#X obj 154 136 *~ 360; +#X obj 63 116 noise~; +#X obj 148 103 phasor~ 0.1; +#X obj 97 194 fft_binaural 128; +#X obj 232 134 nbx 5 14 -1e+37 1e+37 0 0 empty empty empty 0 -8 0 10 +-262144 -1 -1 21 256; +#X msg 24 286 listen_db /yourpath/listen/IRC_1048_C; +#X connect 3 0 6 2; +#X connect 4 0 6 1; +#X connect 5 0 3 0; +#X connect 6 0 2 0; +#X connect 6 1 2 1; +#X connect 7 0 6 3; +#X connect 8 0 6 0; diff --git a/examples/azimuth.gif b/examples/azimuth.gif Binary files differnew file mode 100644 index 0000000..2117182 --- /dev/null +++ b/examples/azimuth.gif diff --git a/examples/coordinates.gif b/examples/coordinates.gif Binary files differnew file mode 100644 index 0000000..7f46de7 --- /dev/null +++ b/examples/coordinates.gif diff --git a/examples/elevation.gif b/examples/elevation.gif Binary files differnew file mode 100644 index 0000000..dd07e0a --- /dev/null +++ b/examples/elevation.gif diff --git a/examples/fft_binaural.pd b/examples/fft_binaural.pd new file mode 100644 index 0000000..72e22ff --- /dev/null +++ b/examples/fft_binaural.pd @@ -0,0 +1,25 @@ +#N canvas 650 270 640 340 10; +#X obj 166 84 inlet~; +#X obj 221 85 inlet~; +#X obj 278 83 inlet~; +#X obj 153 250 outlet~; +#X obj 255 249 outlet~; +#X obj 416 209 block~; +#X obj 402 105 loadbang; +#X msg 413 170 set \$1 1; +#X obj 406 141 expr \$1 / 2; +#X obj 110 83 inlet; +#X obj 151 135 cw_binaural~ \$1 FFT; +#X text 106 22 FFT filtering should be used such than the block size +is equal to half of the considered impulse response length used. Consequently +\, using abstrations similar to this one allows to use spectral filtering +in cw_binaural~ .; +#X connect 0 0 10 0; +#X connect 1 0 10 1; +#X connect 2 0 10 2; +#X connect 6 0 8 0; +#X connect 7 0 5 0; +#X connect 8 0 7 0; +#X connect 9 0 10 0; +#X connect 10 0 3 0; +#X connect 10 1 4 0; diff --git a/matlab_src/README.txt b/matlab_src/README.txt new file mode 100644 index 0000000..a4a7832 --- /dev/null +++ b/matlab_src/README.txt @@ -0,0 +1,3 @@ +This directory contains matlab source code used in our project. +For now, the only available code consists in concerting CIPIC hrtf database, stored in matlab binary file, into wav files that could be used by our system. +The CIPIC hrtf database can be downloaded there: http://interface.idav.ucdavis.edu/data/CIPIC_hrtf_database.zip diff --git a/matlab_src/cipicdb2wav.m b/matlab_src/cipicdb2wav.m new file mode 100755 index 0000000..a5151b1 --- /dev/null +++ b/matlab_src/cipicdb2wav.m @@ -0,0 +1,132 @@ +function cipicdb2wav(src_path, output_path) +%CIPICDB2WAV Summary of this function goes here +% src_path correspond to the root of extracted archive containing +% the CIPIC hrtf database which can be downloaded there +% http://interface.idav.ucdavis.edu/data/CIPIC_hrtf_database.zip +% +% output_path is the directory where the converted database will +% be stored, if the directory does not exists, it will be created + + mkdir(output_path) + + % get a list corresponding the the sets of measures to process + pl = build_processing_list(src_path); + + % those value indexes are given in the cipic documentation, + % they correspond to the azimuths and elevations available in the db + azimuths = [-80 -65 -55 -45:5:45 55 65 80]; + elevations = -45 + 5.625*(0:49); + + % the external uses listen database formalism + % see http://recherche.ircam.fr/equipes/salles/listen/ + % in listen database, angle(90,0) correspond to a point directly to the + % left, and angle(-90,0) correspond to a point directely to the right + % in CIPIC, this is inversed, meaning (90,0) correspond to right + % and so on + % consequenly, we intersed the azimuth angle label to correspond to + % listen formalism + azimuths = -1 * azimuths + + for ihrtf = 1:length(pl) + 'processing' + fname = char(pl(ihrtf,1)) + outdir_name = [output_path '/' char(pl(ihrtf,2))] + mkdir(outdir_name) + data = load(fname) + + for azn=1:length(azimuths) + for eln=1:length(elevations) + + az = azimuths(azn) + el = elevations(eln) + + leftir = squeeze(data.hrir_l(azn, eln,:)); + rightir = squeeze(data.hrir_r(azn, eln,:)); + + stereo(:,1) = leftir; + stereo(:,2) = rightir; + wav_name = sprintf('%s/%s_azim_%.3f_elev_%.3f.wav', outdir_name, char(pl(ihrtf,2)), az, el) + wavwrite(stereo, 44100, 32, wav_name); + end + end + end + +% out_dir='C:\\Users\\david\\Desktop\\\CIPIC_hrtf_database\\rudk_cipic_new'; +% +% azimuths = [-80 -65 -55 -45:5:45 55 65 80]; +% elevations = -45 + 5.625*(0:49); +% +% for ifile = 1:nfiles +% cur_dir = char(dir_list(ifile, :)); +% if cur_dir(1:8) == 'subject_' +% disp(cur_dir) +% cur_dir = strtrim(cur_dir); +% cur_out = [out_dir '\\' cur_dir] +% mkdir(cur_out); +% load([cur_dir '\\\\hrir_final.mat']); +% +% % indexed by azimuth, elevation, and time +% for iazi = 1:length(azimuths) +% for ielev = 1:length(elevations) +% left_data = hrir_l(iazi, ielev, :); +% right_data = hrir_r(iazi, ielev, :); +% [az, el] = normalize_az_el(azimuths(iazi), elevations(ielev)); +% wav_name = sprintf('%s\\%s_az_%d_elev_%.2f.wav', cur_out, cur_dir, az, el) +% stereo_data(:,1) = left_data; +% stereo_data(:,2) = right_data; +% wavwrite(stereo_data, 44100, 32, wav_name); +% if abs(el - 90) < .01 +% wav_name = sprintf('%s\\%s_az_%d_elev_%.2f.wav', cur_out, cur_dir, mod(az+180, 360), el); +% wavwrite(stereo_data, 44100, 32, wav_name); +% end +% end +% end +% end +% end +% +% function [res_az, res_el] = normalize_az_el(az, el) +% % elevation should be between - 90 and + 90 +% if el > 226 %% HACK +% +% 'hack' +% +% el = -(el - 225); +% az = az + 180; +% +% end +% +% if el > 90 + .005 +% el = 180 -el; +% az = 180 + az; +% end +% +% +% az = mod(az, 360); +% if (az < 0) +% az = az + 360; +% end +% +% res_az = az; +% res_el = el; + function res = build_processing_list(path) + res = [] + c = {'tutu','tata','titi'} + d = [c; {'tutu2','tata2','titi2'}] + e = [d; {'tutu3','tata4','titi5'}] + path + l = dir([path '/standard_hrir_database/subject_*']) + for i=1:length(l) + fullname = [path '/standard_hrir_database/' l(i).name '/hrir_final.mat'] + res = [res; {fullname, char(l(i).name)}] + end +% those files do not contains the same fields +% p = [path '/special_kemar_hrir/kemar_frontal/large_pinna_frontal.mat'] +% res = [res; {p, 'kemar_frontal_large_pinna'}] +% p = [path '/special_kemar_hrir/kemar_frontal/small_pinna_frontal.mat'] +% res = [res; {p, 'kemar_frontal_small_pinna'}] +% p = [path '/special_kemar_hrir/kemar_horizontal/large_pinna_final.mat'] +% res = [res; {p, 'kemar_horizontal_large_pinna'}] +% p = [path '/special_kemar_hrir/kemar_horizontal/small_pinna_final.mat'] +% res = [res; {p, 'kemar_horizontal_small_pinna'}] + end +end
\ No newline at end of file diff --git a/src/binaural_processor.cpp b/src/binaural_processor.cpp new file mode 100644 index 0000000..22a8597 --- /dev/null +++ b/src/binaural_processor.cpp @@ -0,0 +1,164 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + +#include "binaural_processor.hpp" +#include "riff.hpp" +#include "fft_riff.hpp" +#include "logstring.hpp" +#include "flyweight_ir_factory.hpp" + + + +BinauralProcessor::~BinauralProcessor() +{ + delete _leftdelay; + delete _rightdelay; + delete _leftriff; + delete _rightriff; +} + +int BinauralProcessor::set_listen_db(string path) +{ + string re("IRC_[0-9]+_C_R[0-9]+_T([-+]?[0-9]*\\.?[0-9]+)_P([-+]?[0-9]*\\.?[0-9]+)\\.wav"); + return set_hrtf(path, re, true, true); +} + +int BinauralProcessor::set_cipic_db(string path) +{ + string re("subject_[0-9]+_azim_([-+]?[0-9]*\\.?[0-9]+)_elev_([-+]?[0-9]*\\.?[0-9]+)\\.wav"); + return set_hrtf(path, re, true, false); +} + + +int BinauralProcessor::set_hrtf(string hrtf_path, string regex, bool azfirst, bool vertpolar) +{ + ir_key k; + + // slog << "Setting hrtf " << hrtf_path << endl; + k.path = hrtf_path; + k.length = _length; + k.iir_regex = regex; + k.is_azimuth_first = azfirst; + k.spectral = _spectral_processing; + k.minp_ap_dec = _minphase_allpass_decomposition; + k.vertical_polar_coords = vertpolar; + + _newhrtfcont = FlyweightIrFactory::instance()->hrtf_set_get(k); + _itdcont = FlyweightIrFactory::instance()->itd_set_get(k); + + // slog << "hrtfcont" << _newhrtfcont << endl; + // slog << "itdcont" << _itdcont << endl; + + return 0; +} + + + +void BinauralProcessor::process(float* input, float* azimuths, float* elevations, float* left_output, float* right_output, int n) +{ + // FIXME: BUG IF LEFT OUTPUT = INPUT.... + if (!_newhrtfcont) + { + // if no hrtf is selected, output a null signal + for (int i = 0; i < n; ++i) + left_output[i] = right_output[i] = 0; + return; + } + + // for a given block, consider only the last requested azimuth and elevation + float az = azimuths[n-1]; + float el = elevations[n-1]; + + // update ic to store identifiers of impulse responses to be interpolated + interp_cdts ic; + _newhrtfcont->set_candidates(ic, az, el); + + _newhrtfcont->update_from_candidates(ic, _leftriff->coeff_get(), _rightriff->coeff_get()); + + // Apply Fractional Delay + if (_itdcont) + { + float itd = _itdcont->itd_from_candidates(ic); + _leftdelay->process(input, left_output, itd < 0 ? -itd: 0, n); + _rightdelay->process(input, right_output, itd > 0 ? itd : 0, n); + } + else + for (int i = 0; i < n; ++i) + left_output[i] = right_output[i] = input[i]; + + // Filter the signal + _leftriff->process(left_output, left_output, n); + _rightriff->process(right_output, right_output, n); +} + + + + +BinauralProcessor::BinauralProcessor(int impulse_response_size, string& filtering_method, string& delay_method): _newhrtfcont(NULL), _itdcont(NULL) +{ + + //slog << "New binaural processor using IR size " << impulse_response_size << ", filter method " << filtering_method << ", delay method " << delay_method << endl; + + // check if the requested filtering method is supported + if (!filtering_method.compare("RIFF")) + _spectral_processing = false; + else if (!filtering_method.compare("FFT")) + _spectral_processing = true; + else + { + slog << filtering_method << "is not a valid filtering method name: try \"FFT\" or \"RIFF\"" << endl; + throw 0; + } + + // check if the requested impulse response size is supported + // FIXME: check if it is a correct test!! => PUT THAT IN THE FLYWEIGHT + _length = impulse_response_size; + while (impulse_response_size > 2) + impulse_response_size /= 2; + if (impulse_response_size != 2) + { + slog << "Impulse response size should be a positive power of 2: recommended values are 128, 256 or 512" << endl; + throw 0; + } + + // ********************* + // TODO + // at init, use only a dirac impulsion?? + + if (_spectral_processing) + { + _leftriff = new FftRiff(_length); + _rightriff = new FftRiff(_length); + } + else + { + _leftriff = new Riff(_length); + _rightriff = new Riff(_length); + } + + // delays + _leftdelay = Delay::create(delay_method, _length); + _rightdelay = Delay::create(delay_method, _length); + _minphase_allpass_decomposition = delay_method.compare("nodelay"); +} diff --git a/src/binaural_processor.hpp b/src/binaural_processor.hpp new file mode 100644 index 0000000..8106797 --- /dev/null +++ b/src/binaural_processor.hpp @@ -0,0 +1,73 @@ +/*
+ cw_binaural~: a binaural synthesis external for pure data
+ by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan
+ and Anne Sedes - sedes.anne@gmail.com
+ Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 +
+
+ This program is free software: you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+*/
+
+
+#ifndef BINAURALPROCESSOR_HPP_
+# define BINAURALPROCESSOR_HPP_
+
+# include <string>
+# include "generic_riff.hpp"
+# include "delay.hpp"
+# include "ir_key.hpp"
+# include "hrtfcont.hpp"
+# include "itdcont.hpp"
+
+
+class BinauralProcessor
+{
+public:
+ BinauralProcessor(int impulse_response_size, std::string& filtering_method, std::string& delay_method);
+ ~BinauralProcessor();
+
+ // do the processing for a given filter
+ void process(float* input, float *azimuths, float* elevations, float* left_output, float* right_output, int n);
+
+ // load a HRTF database, assuming it is LISTEN database
+ int set_listen_db(string path);
+
+ // load a HRTF database, assuming it is CIPIC database
+ int set_cipic_db(string path);
+
+ // load any hrtf database
+ int set_hrtf(string hrtf_path, string regex, bool azfirst, bool vertpolar);
+
+protected:
+
+
+ // we should do an optimized riff structure
+ GenericRiff *_leftriff, *_rightriff;
+ Delay* _leftdelay, *_rightdelay;
+
+ // creation arguments
+ int _length;
+ bool _spectral_processing;
+ bool _minphase_allpass_decomposition;
+
+ // hrtf ad itd database
+ HrtfCont *_newhrtfcont;
+ ItdCont *_itdcont;
+};
+
+
+#endif
diff --git a/src/cw_binaural~.cpp b/src/cw_binaural~.cpp new file mode 100644 index 0000000..03468f0 --- /dev/null +++ b/src/cw_binaural~.cpp @@ -0,0 +1,237 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + +#include <m_pd.h> +#include "binaural_processor.hpp" +#include "logstring.hpp" + +static t_class *cw_binaural_tilde_class; + +typedef struct _cw_binaural_tilde +{ + t_object x_obj; + BinauralProcessor *bp; + float default_input; +} t_cw_binaural_tilde; + + +t_int *cw_binaural_tilde_perform(t_int *w) +{ + t_cw_binaural_tilde *x = (t_cw_binaural_tilde *)(w[1]); + + // 3 signal input: signal, azimuth and elevation + t_sample *input = (t_sample *)(w[2]); + t_sample *azimuths = (t_sample *)(w[3]); + t_sample *elevations = (t_sample *)(w[4]); + + // 2 signal outputs: for the left and right ear + t_sample *left_output = (t_sample *)(w[5]); + t_sample *right_output = (t_sample *)(w[6]); + + int n = (int)(w[7]); + + x->bp->process(input, azimuths, elevations, left_output, right_output, n); + + return (w+8); +} + + + + + +void cw_binaural_tilde_dsp(t_cw_binaural_tilde *x, t_signal **sp) +{ + dsp_add(cw_binaural_tilde_perform, 7, x, sp[0]->s_vec, sp[1]->s_vec, + sp[2]->s_vec, sp[3]->s_vec, sp[4]->s_vec, + sp[0]->s_n); +} + +// called when the message listen_db is sent to the object +void cw_binaural_set_listen_db(t_cw_binaural_tilde *x, t_symbol *hrtf_path) +{ + try + { + x->bp->set_listen_db(string(hrtf_path->s_name)); + } + catch (int n) {} + logstring2pdterm(); +} + +// called when the message cipic_db is sent to the object +void cw_binaural_set_cipic_db(t_cw_binaural_tilde *x, t_symbol *hrtf_path) +{ + try + { + x->bp->set_cipic_db(string(hrtf_path->s_name)); + } + catch (int n) {} + logstring2pdterm(); +} + +// convert a symbol to boolean +static bool symb2bool(t_symbol* symb) +{ + string s(symb->s_name); + + for (size_t i = 0; i < s.length(); ++i) + s[i] = toupper(s[i]); + return s == "TRUE"; +} + +// called when the message sethrtf is sent to the object +// message allowing to load any hrtf database +// hrtf_path correspond to the path to the directory storing the hrtf as wav files +// fname_regex is a regex representing the name of files containing impulse responses +// it should contain a group for each parameter to extract (azimuth and elevation) +// az_first tells wether the first group corresponds to azimuth or not +// vertical_polar_coords is set to true if the database is expressed in vertical polar +// coordinates, if set to false, it refers interaural polar coordinates +void cw_binaural_set_hrtf_db(t_cw_binaural_tilde *x, + t_symbol *symb_hrtf_path, + t_symbol *symb_fname_regex, + t_symbol *symb_az_first, + t_symbol *symb_vertical_polar_coords) +{ + string hrtf_path(symb_hrtf_path->s_name); + string fname_regex(symb_fname_regex->s_name); + bool az_first = symb2bool(symb_az_first); + bool vertical_polar_coords = symb2bool(symb_vertical_polar_coords); + + // Pd does not allows to transmit backslashes in messages, + // so we replace each occurence of double slashes with a back slash + size_t found; + while ((found = fname_regex.find("//")) != string::npos) + fname_regex.replace(found, 2, 1, '\\'); + + try + { + x->bp->set_hrtf(hrtf_path, fname_regex, az_first, vertical_polar_coords); + } + catch (int n) {} + logstring2pdterm(); +} + + + +void cw_binaural_tilde_free(t_cw_binaural_tilde *x) +{ + if (x->bp) + delete x->bp; + logstring2pdterm(); +} + +void *cw_binaural_tilde_new(t_symbol *obj_name, int argc, t_atom *argv) + +{ + // default creation arguments + int impulse_response_size = 128; + std::string filtering_method = std::string("RIFF"); + std::string delay_method = std::string("Hermite4"); + + // parsing creation arguments + switch (argc < 3 ? argc : 3) + { + case 3: + delay_method = std::string(atom_getsymbol(argv + 2)->s_name); + case 2: + filtering_method = std::string(atom_getsymbol(argv + 1)->s_name); + case 1: + impulse_response_size = (int) atom_getfloat(argv); + default: + break; + } + + // creating object + t_cw_binaural_tilde *x = (t_cw_binaural_tilde *)pd_new(cw_binaural_tilde_class); + x->default_input = 0; + x->bp = NULL; + + try + { + // call to the constructor + x->bp = new BinauralProcessor(impulse_response_size, filtering_method, delay_method); + } + catch (int n) + { + x->bp = 0; + logstring2pdterm(); + return (void*) 0; + } + + // azimuth inlet + inlet_new(&x->x_obj, &x->x_obj.ob_pd, &s_signal, &s_signal); + // elevation inlet + inlet_new(&x->x_obj , &x->x_obj.ob_pd, &s_signal, &s_signal); + + // left outlet + outlet_new(&x->x_obj, &s_signal); + // right outlet + outlet_new(&x->x_obj, &s_signal); + + logstring2pdterm(); + + return (void *)x; +} + + +extern "C" +#ifdef WIN32 +__declspec(dllexport) void cw_binaural_tilde_setup(void) +#else + void cw_binaural_tilde_setup(void) +#endif +{ + // creation of the cw_binaural~ instance + cw_binaural_tilde_class = + class_new(gensym("cw_binaural~"), + (t_newmethod)cw_binaural_tilde_new, + (t_method) cw_binaural_tilde_free, sizeof(t_cw_binaural_tilde), + CLASS_DEFAULT, + A_GIMME, 0); + + // sound processing + class_addmethod(cw_binaural_tilde_class, + (t_method)cw_binaural_tilde_dsp, gensym("dsp"), (t_atomtype) 0); + + // management of messages used to set a listen database + class_addmethod(cw_binaural_tilde_class, + (t_method)cw_binaural_set_listen_db, + gensym("listen_db"), + A_DEFSYMBOL, (t_atomtype) 0); + + // management of messages used to set a CIPIC database + class_addmethod(cw_binaural_tilde_class, + (t_method)cw_binaural_set_cipic_db, + gensym("cipic_db"), + A_DEFSYMBOL, (t_atomtype) 0); + + // management of messages used to set any database + class_addmethod(cw_binaural_tilde_class, + (t_method)cw_binaural_set_hrtf_db, + gensym("set_hrtf_db"), + A_DEFSYMBOL, A_DEFSYMBOL, A_DEFSYMBOL, A_DEFSYMBOL, (t_atomtype) 0); + + CLASS_MAINSIGNALIN(cw_binaural_tilde_class, t_cw_binaural_tilde, default_input); +} + diff --git a/src/delay.cpp b/src/delay.cpp new file mode 100644 index 0000000..047f043 --- /dev/null +++ b/src/delay.cpp @@ -0,0 +1,162 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + +#include <math.h> +#include "delay.hpp" +#include "logstring.hpp" + +Delay* Delay::create(string delay_name, unsigned size) +{ + //slog << "creating delay " << delay_name << endl; + // kind of factory design pattern + if (!delay_name.compare("Hermite4")) + return new SCDelay(size); + else if (!delay_name.compare("nodelay")) + return new NoDelay(size); + else if (!delay_name.compare("linear")) + return new BasicInterpolDelay(size); + else if (!delay_name.compare("nofractional")) + return new DummyDelay(size); + else if (!delay_name.compare("6points")) + return new CZDelay(size); + else + { + slog << delay_name << " delay method is not supported" << endl; + throw 0; + } +} + +Delay::Delay(unsigned size): _size(size) +{ + // FIXME + // ASUMPTION: ITD WILL NEVER BE > SIZE/2 + // IT SHOULD BE CHECKED DURING ITD COMPUTATION + _fifo = new float[size]; + _fifo_id = 0; + for (unsigned i = 0; i < size; ++i) + _fifo[i] = 0; + _cur_delay = 0; +} + +Delay::~Delay() +{ + delete [] _fifo; +} + +inline void Delay::push_value(float f) +{ + _fifo_id = (_fifo_id ? _fifo_id - 1 : _size - 1); + _fifo[_fifo_id] = f; +} + +// That ugly macro has been done because of the inneficiency +// of virtual interpolation functions for demanding RT purposes +#define DELAY_PROCESS(code) \ + unsigned i; \ + _old_delay = _cur_delay; _cur_delay = new_delay; \ + float delay = _old_delay; \ + const float inc_delay = (_cur_delay - _old_delay) / n; \ + \ + for (i = 0; i < n; ++i, delay += inc_delay) \ + { \ + push_value(input[i]); \ + /* compute the id of the current output element \ + delayed by and integer value */ \ + float floor_delay = floorf(delay); \ + \ + /* get the current element id in the fifo */ \ + unsigned id = _fifo_id + (unsigned) floor_delay; \ + if (id >= _size) \ + id -= _size; \ + \ + /* compute the fractional part of the delay */ \ + float frac_del = delay - floor_delay; \ + \ + code \ + \ + output[i] = res; \ + } + + +void DummyDelay::process(const float* input, float* output, float new_delay, unsigned n) +{ + DELAY_PROCESS( + float res = _fifo[id]; + ) +} + +void BasicInterpolDelay::process(const float* input, float* output, float new_delay, unsigned n) +{ + DELAY_PROCESS( + unsigned older = (id + 1 == _size ? 0: id + 1); + float res = _fifo[id] * frac_del + (1 - frac_del) * _fifo[older]; + ) +} + +void SCDelay::process(const float* input, float* output, float new_delay, unsigned n) +{ + DELAY_PROCESS( + float y0 = _fifo[id]; + float y1 = _fifo[(id + 1) % _size]; + float y2 = _fifo[(id + 2) % _size]; + float y3 = _fifo[(id + 3) % _size]; + + // 4-point, 3rd-order Hermite (x-form) + float c0 = y1; + float c1 = 0.5 * (y2 - y0); + float c2 = y0 - 2.5 * y1 + 2. * y2 - 0.5 * y3; + float c3 = 0.5 * (y3 - y0) + 1.5 * (y1 - y2); + + float res = ((c3 * frac_del + c2) * frac_del + c1) * frac_del + c0; + ) +} + +void CZDelay::process(const float* input, float* output, float new_delay, unsigned n) +{ + // 6 points interpolation + DELAY_PROCESS( + float ym2 = _fifo[id]; + float ym1 = _fifo[(id + 1) % _size]; + float y0 = _fifo[(id + 2) % _size]; + float y1 = _fifo[(id + 3) % _size]; + float y2 = _fifo[(id + 4) % _size]; + float y3 = _fifo[(id + 5) % _size]; + + float a0= y0; + float a1= (1./12.)*ym2 - (2./3.)*ym1 + (2./3.)*y1 - (1./12.)*y2; + float a2= (-1./24.)*ym2 + (2./3.)*ym1 - (5./4.)*y0 + (2./3.)*y1 - (1./24.)*y2; + float a3= (-3./8.)*ym2 + (13./8.)*ym1 - (35./12.)*y0 + (11./4.)*y1 - (11./8.)*y2 + (7./24.)*y3; + float a4= (13./24.)*ym2 - (8./3.)*ym1 + (21./4.)*y0 - (31./6.)*y1 + (61./24.)*y2 - (1./2.)*y3; + float a5= (-5./24.)*ym2 + (25./24.)*ym1 - (25./12.)*y0 + (25./12.)*y1 - (25./24.)*y2 + (5./24.)*y3; + + float res = a0 + frac_del * (a1 + frac_del * (a2 + frac_del * (a3 + frac_del * (a4 + frac_del * a5)))); + ) +} + + +void NoDelay::process(const float *input, float *output, float new_delay, unsigned n) +{ + for (unsigned i = 0; i < n; ++i) + output[i] = input[i]; +} diff --git a/src/delay.hpp b/src/delay.hpp new file mode 100644 index 0000000..c6c0a05 --- /dev/null +++ b/src/delay.hpp @@ -0,0 +1,91 @@ +/*
+ cw_binaural~: a binaural synthesis external for pure data
+ by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan
+ and Anne Sedes - sedes.anne@gmail.com
+ Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 +
+
+ This program is free software: you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+*/
+
+
+#ifndef DELAY_HPP_
+# define DELAY_HPP_
+
+# include <string>
+
+using namespace std;
+
+class Delay
+{
+public:
+ Delay(unsigned size);
+ virtual ~Delay();
+ virtual void process(const float* input, float* output, float new_delay, unsigned n) = 0;
+ static Delay* create(string delay_name, unsigned size);
+protected:
+ inline void push_value(float f);
+
+ const unsigned _size;
+ unsigned _fifo_id;
+ float * _fifo;
+ float _old_delay, _cur_delay;
+};
+
+
+// Does no fractional Delay
+class DummyDelay: public Delay
+{
+public:
+ DummyDelay(int size): Delay(size) {}
+ virtual void process(const float* input, float* output, float new_delay, unsigned n);
+};
+
+// 1st order fractional Delay
+class BasicInterpolDelay: public Delay
+{
+public:
+ BasicInterpolDelay(int size): Delay(size) {}
+ virtual void process(const float* input, float* output, float new_delay, unsigned n);
+};
+
+// SuperCollider fractional Delay
+class SCDelay: public Delay
+{
+public:
+ SCDelay(int size): Delay(size) {}
+ virtual void process(const float* input, float* output, float new_delay, unsigned n);
+};
+
+
+// 6th order fractional Delay
+class CZDelay: public Delay
+{
+public:
+ CZDelay(int size): Delay(size) {}
+ virtual void process(const float* input, float* output, float new_delay, unsigned n);
+};
+
+
+class NoDelay: public Delay
+{
+public:
+ NoDelay(int size): Delay(size) {}
+ virtual void process(const float* input, float* output, float new_delay, unsigned n);
+};
+
+#endif
diff --git a/src/fft_riff.cpp b/src/fft_riff.cpp new file mode 100644 index 0000000..bce1fd1 --- /dev/null +++ b/src/fft_riff.cpp @@ -0,0 +1,90 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include <math.h> +#include <m_pd.h> + +#include "fft_riff.hpp" +#include "logstring.hpp" + + +FftRiff::FftRiff(int size): GenericRiff(size) +{ + delete [] _coeffs; + _coeffs = new float[2* size]; + + for (int i = 0; i < 2*size; ++i) + _coeffs[i] = 0; + + _tmp_input = new float[2 * size]; + _tmp_output = new float[2 * size]; + _tmp_rest = new float[2*size]; + // we divide by size for not having to normalize after a fft + ifft + for (int i = 0; i < 2*size; ++i) + _tmp_rest[i] = 0; +} + +FftRiff::~FftRiff() +{ + delete [] _tmp_input; + delete [] _tmp_output; + delete [] _tmp_rest; +} + + +void FftRiff::process(float* input, float* output, int n) +{ + int i; + + for (i = 0; i < n; ++i) + _tmp_input[i] = input[i] / (_size * 2); + for (i = n; i < _size * 2; ++i) + _tmp_input[i] = 0; + + mayer_realfft(2* _size, _tmp_input); + + _tmp_output[0] = _tmp_input[0] * _coeffs[0]; + _tmp_output[_size] = _tmp_input[_size] * _coeffs[_size]; + + /// multiply the spectrum of the input and the hrtf + for (int i = 1; i < _size; ++i) + { + // real part + _tmp_output[i] = _tmp_input[i] * _coeffs[i] - _tmp_input[2*_size-i] * _coeffs[2*_size-i]; + // imag part + _tmp_output[2*_size-i] = _tmp_input[i] * _coeffs[2*_size-i] + _tmp_input[2*_size-i] * _coeffs[i]; + } + + + mayer_realifft(2*_size, _tmp_output); + + // output the current block beginning + the corresponding rests + for (i = 0; i < n; ++i) + output[i] = _tmp_output[i] + _tmp_rest[i]; + for (i = n; i < 2*_size-n; ++i) + _tmp_rest[i-n] = _tmp_rest[i] + _tmp_output[i]; + for (i = 2*_size - n; i < 2 * _size; ++i) + _tmp_rest[i-n] = _tmp_output[i]; +} diff --git a/src/fft_riff.hpp b/src/fft_riff.hpp new file mode 100644 index 0000000..ce658b5 --- /dev/null +++ b/src/fft_riff.hpp @@ -0,0 +1,54 @@ +/*
+ cw_binaural~: a binaural synthesis external for pure data
+ by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan
+ and Anne Sedes - sedes.anne@gmail.com
+ Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 +
+
+ This program is free software: you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+*/
+
+
+#ifndef FFTRIFF_HPP_
+# define FFTRIFF_HPP_
+
+# include "generic_riff.hpp"
+
+class FftRiff: public GenericRiff
+{
+public:
+ FftRiff(int size);//: GenericRiff(size) {}
+ ~FftRiff();
+ virtual void process(float* input, float* output, int n);
+ //float* coeff_get() {return _coeffs;}
+protected:
+ // void process_core(float* input, float*output);
+ //void swap_fpointer((float*)& a, (float*)& b);
+ //const int _size;
+ //int _fifo_id;
+ //float *_last_output;
+ //float *_window;
+ float *_tmp_input;
+ float *_tmp_output;
+ float *_tmp_rest;
+
+ //float *_cur_input, *_next_input;
+ //float *_cur_output, *_next_output;
+ //int _inputpos;
+};
+
+#endif
diff --git a/src/flyweight_ir_factory.cpp b/src/flyweight_ir_factory.cpp new file mode 100644 index 0000000..df29480 --- /dev/null +++ b/src/flyweight_ir_factory.cpp @@ -0,0 +1,121 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include "flyweight_ir_factory.hpp" +#include "logstring.hpp" +#include "raw_wav_hrtfcont.hpp" +#include "unprocessed_fixed_size_hrtfcont.hpp" +#include "spectral_hrtfcont.hpp" +#include "minphase_hrtfcont.hpp" + +FlyweightIrFactory* FlyweightIrFactory::_instance = NULL; + +FlyweightIrFactory* FlyweightIrFactory::instance() +{ + //slog << "flyw instance" << endl; + if (!_instance) + _instance = new FlyweightIrFactory(); + return _instance; +} + +HrtfCont* FlyweightIrFactory::hrtf_set_get(const ir_key& k) +{ + //slog << "calling the factory" << endl; + // if the key is in the map, then it is valid, and an be returned as it + if (_hrtf_db.find(k) != _hrtf_db.end()) + return _hrtf_db[k]; + + // slog << "key not in database, creating new db" << endl; + + // check that the key is valid + if (!is_hrtf_set_key_consistent(k)) + { + slog << "key not consistent" << endl; + return NULL; + } + + // By convention, length == 0 means + // getting the unprocessed hrtf at its original size + if (!k.length) + return (_hrtf_db[k] = new RawWavHrtfCont(k)); + + // spectral representation of the filter database + // to be used in overlap-add + if (k.spectral) + return (_hrtf_db[k] = new SpectralHrtfCont(k)); + + // minmum phase filter + if (k.minp_ap_dec) + return (_hrtf_db[k] = new MinPhaseHrtfCont(k)); + + // Unprocessed impulse response of fixed size + // this test is useless + if (!k.minp_ap_dec && !k.spectral) + return (_hrtf_db[k] = new UnprocessedFixedSizeHrtfCont(k)); + + //slog << "returb null" << endl; + return NULL; +} + +bool FlyweightIrFactory::is_hrtf_set_key_consistent(const ir_key& k) +{ + if (!k.path.length() || !k.iir_regex.length()) + return false; + + if (!k.length && (k.minp_ap_dec || k.spectral)) + return false; + + if ((k.minp_ap_dec || k.spectral) && !ispowerof2(k.length)) + { + slog << "For spectral filtering or for minphase decompostion, " + << "the impulse response size should be a power of 2" << endl; + return false; + } + return true; +} + +bool FlyweightIrFactory::ispowerof2(size_t n) +{ + return (n > 0) && !(n & (n - 1)); +} + +ItdCont* FlyweightIrFactory::itd_set_get(ir_key k) +{ + // if no minphase/allpass decomposition has been done + // then, we won't use fractional delay and ITD information + if (!k.minp_ap_dec) + return NULL; + + // the itd will be the same idependently of the spectral field + // so we put spectral to an arbitrary value to avoid duplicates + k.spectral = false; + + // if the key is in the map, then it is valid, and an be returned as it + if (_itd_db.find(k) != _itd_db.end()) + return _itd_db[k]; + + // Only one ITD computation method is available for now + return (_itd_db[k] = new ItdCont(k)); +} diff --git a/src/flyweight_ir_factory.hpp b/src/flyweight_ir_factory.hpp new file mode 100644 index 0000000..3a712eb --- /dev/null +++ b/src/flyweight_ir_factory.hpp @@ -0,0 +1,71 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#ifndef FLYWEIGHT_IR_FACTORY_HPP_ +# define FLYWEIGHT_IR_FACTORY_HPP_ + +# include <map> +# include "ir_key.hpp" + +# include "hrtfcont.hpp" +# include "itdcont.hpp" + + +// This class is a mixed flyweight and factory approach +// aimed to obtain the requested impulse response +class FlyweightIrFactory +{ +public: + // there should be only one flyweightfactory + // so we use a singleton design pattern + static FlyweightIrFactory* instance(); + + // get a hrtf set container corresponding + // to a given set of parameters + HrtfCont* hrtf_set_get(const ir_key& k); + // get a set of itd corresponding to a given set of parameters + ItdCont* itd_set_get(ir_key k); +private: + // check if a given key is consistent + bool is_hrtf_set_key_consistent(const ir_key& k); + + // tells if a number is a power of 2 + bool ispowerof2(size_t n); + + // All Hrtf sets will be stored in a hash map + // indexed by the requested preprocessing options + // on the hrtf set + map<ir_key, HrtfCont*, ir_key_isless> _hrtf_db; + + // a separate map is used for ITD computations + // on a longer term, it will allow to compute itd using + // different methods + map<ir_key, ItdCont*, ir_key_isless> _itd_db; + + FlyweightIrFactory() {} + static FlyweightIrFactory* _instance; +}; + +#endif diff --git a/src/generic_riff.cpp b/src/generic_riff.cpp new file mode 100644 index 0000000..619e6ed --- /dev/null +++ b/src/generic_riff.cpp @@ -0,0 +1,38 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + +#include "generic_riff.hpp" + +GenericRiff::GenericRiff(int size): _size(size) +{ + _coeffs = new float[size]; + + for (int i = 0; i < size; ++i) + _coeffs[i] = 0; +} + +GenericRiff::~GenericRiff() +{ + delete [] _coeffs; +} diff --git a/src/generic_riff.hpp b/src/generic_riff.hpp new file mode 100644 index 0000000..c252011 --- /dev/null +++ b/src/generic_riff.hpp @@ -0,0 +1,45 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + +#ifndef GENERICRIFF_HPP_ +# define GENERICRIFF_HPP_ + +// this is an astract class for Riff filters +// TODO: REMOVE the size attribute: it may change +// depending on the block size +// currently is is dependant of the HRTF impulse response size... + +class GenericRiff +{ +public: + GenericRiff(int size); + virtual ~GenericRiff(); + virtual void process(float* input, float* output, int n) = 0; + float* coeff_get() {return _coeffs;} +protected: + const int _size; + float *_coeffs; +}; + +#endif diff --git a/src/hrtfcont.cpp b/src/hrtfcont.cpp new file mode 100644 index 0000000..7e018f0 --- /dev/null +++ b/src/hrtfcont.cpp @@ -0,0 +1,256 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include <iostream> +#include <math.h> +#include "hrtfcont.hpp" +#include "logstring.hpp" + +#define DEG2RAD (M_PIl/180) +#define RAD2DEG (180/M_PIl) +#define EPS 0.0000000000001 + +HrtfCont::HrtfCont(const ir_key& k): + _vert_pol_coords(k.vertical_polar_coords) {} + +// set elevation in [-90, +90], and azimuth in [-180, 180[ +void HrtfCont::normalize_vertpolar_coords(float& azimuth, float& elevation) const +{ + // convert any elevation into a usefull value + // the value used should lie between -90 and +90 + elevation = fmod(elevation, 360); + if (elevation < 0) + elevation += 360; + + if (elevation > 270) + elevation -= 360; + else if (elevation > 90) + { + elevation = 180 - elevation; + azimuth += 180; + } + + // azimuth values should be in [0,360[ + azimuth = fmod(azimuth, 360); + if (azimuth < -180) + azimuth += 360; + else if (azimuth >= 180) + azimuth -= 360; +} + + +// return the distance in degree between 2 normalized angles +// ie both angles are in [0, 360[ +float HrtfCont::angular_dist(float a1, float a2) const +{ + float d = fabs(a1-a2); + return d <= 180 ? d : 360 - d; +} + + +// iangle1, iangle2 +void HrtfCont::add_a2_candidates(interp_cdts& ic, float a1_key, float a2, float weight) +{ + // get impulse responses corresponding to the current elevation + const angle2_cont& a2_candidates = _m[a1_key]; + //slog << "interp for fixed angle1=" << a1_key << ", optimal angle2="<< a2<< ", weight=" << weight << endl; + // there is only one angle2 measure at index angle1 + if (a2_candidates.size() == 1) + return ic.add(a1_key, a2_candidates.begin()->first, weight); + + // iterator the the 1st element of key bigger or equal + angle2_cit supeq_a2it = a2_candidates.lower_bound(a2); + + // the requested angle2 corresponds exactly to an available measure + if (supeq_a2it != a2_candidates.end() && a2 == supeq_a2it->first) + return ic.add(a1_key, supeq_a2it->first, weight); + + // index candidates for angle2 + float a2_cand1, a2_cand2; + if (supeq_a2it == a2_candidates.end() || supeq_a2it == a2_candidates.begin()) + { + //slog << "angle2 btwn 2 extremes" << (supeq_a2it == a2_candidates.end()) << (supeq_a2it == a2_candidates.begin()) << endl; + // requested angle2 is bigger than all available measures + // or smaller than all available measures + // interpolation between extreme measures can be done + // since angle2 is in [-180, 180[ + a2_cand1 = a2_candidates.begin()->first; + supeq_a2it = a2_candidates.end(); + advance(supeq_a2it, -1); + a2_cand2 = supeq_a2it->first; + } + else + { + // general case + //slog << "angle2 interp general case" << endl; + a2_cand1 = supeq_a2it->first; + advance(supeq_a2it, -1); + a2_cand2 = supeq_a2it->first; + } + // angular distance between the 2 candidates + float a2_cands_dist = angular_dist(a2_cand1, a2_cand2); + // weight of the second candidate + float a2_cand2_weight = angular_dist(a2, a2_cand1) / a2_cands_dist; + //slog << "dist between a2 candidates" << a2_cands_dist << endl; + //slog << "angular_dist(a2, a2_cand1)" << angular_dist(a2, a2_cand1) << endl; + + // add 2 candidates + ic.add(a1_key, a2_cand1, weight * (1-a2_cand2_weight)); + ic.add(a1_key, a2_cand2, weight * a2_cand2_weight); +} + + +void HrtfCont::set_candidates(interp_cdts& ic, float az, float el) +{ + // assumption: ic.size == 0, and map size != 0, not checked for RT issues :-( + + // index angles in the hrtf map + float iangle1, iangle2; + + //slog << "set candidates az:" << az << ", el:" << el << endl; + // get the angular indexes expressed in the coordianates of the HRTF db to use + if (_vert_pol_coords) + { + // Database using vertical polar coordinates + normalize_vertpolar_coords(az, el); + iangle1 = el; + iangle2 = az; + } + else + { + // Database using interaural polar coordiantes + vertpol2interaurpol(az,el); + iangle1 = az; + iangle2 = el; + } + //slog << "indexes " << iangle1 << ", " << iangle2 << endl; + + // there is only one available value for first index + if (_m.size() == 1) + return add_a2_candidates(ic, _m.begin()->first, iangle2, 1); + + // get the first index value >= iangle1 + angle1_cit supeq_a1it = _m.lower_bound(iangle1); + + // the requested angle1 corresponds exactly to an indexed element + if (supeq_a1it != _m.end() && iangle1 == supeq_a1it->first) + return add_a2_candidates(ic, iangle1, iangle2, 1); + + if (supeq_a1it == _m.begin() || supeq_a1it == _m.end()) + { + // requested angle1 is strictly smaller than + // the smallest index of the database + // or strictly bigger than the biggest index of the db + float range_without_measures; + float w2; + if (supeq_a1it == _m.end()) + { + // requested index1 is bigger than available indexes + advance(supeq_a1it, -1); + range_without_measures = 2*(90 - supeq_a1it->first); + w2 = iangle1 - supeq_a1it->first; + } + else + { + // requested index1 is smaller than available indexes + range_without_measures = -2*(-90 - supeq_a1it->first); + w2 = supeq_a1it->first - iangle1; + } + w2 /= range_without_measures; + add_a2_candidates(ic, supeq_a1it->first, iangle2, 1-w2); + iangle2 = iangle2 > 0 ? iangle2 - 180 : iangle2 + 180; + add_a2_candidates(ic, supeq_a1it->first, iangle2, w2); + } + else + { + // general case: intepolation using 2 different indexed elements + float index1_candidate1, index1_candidate2; + index1_candidate2 = supeq_a1it->first; + advance(supeq_a1it, -1); + index1_candidate1 = supeq_a1it->first; + // distance between index1 candidates. NB: candidate 2 is bigger than candidate 1! + float dcands = index1_candidate2 - index1_candidate1; + add_a2_candidates(ic, index1_candidate2, iangle2, (iangle1-index1_candidate1)/dcands); + add_a2_candidates(ic, index1_candidate1, iangle2, (index1_candidate2-iangle1)/dcands); + } +} + + +// TODO: could be optimized using templates +void HrtfCont::update_from_candidates(const interp_cdts& ic, float* left, float* right) +{ + //slog << "update from candidates " << left << right << endl; + const float *lc, *rc; + float wc; + ir_buffer* irb; + + irb = &(_m[ic.angle_index1[0]][ic.angle_index2[0]]); + //slog << ic.angle_index1[0] << " " << ic.angle_index2[0] << " " << irb->fname << endl; + // return; + + wc = ic.weight[0]; + lc = irb->lbuf; + rc = irb->rbuf; + //slog << "update from " << _m[ic.angle_index1[0]][ic.angle_index2[0]].fname << " " << wc << endl; + //slog << lc << rc << left << right << endl; + + for (size_t i = 0; i < _ir_length; ++i) + { + left[i] = wc * lc[i]; + right[i] = wc * rc[i]; + } + + for (size_t icand = 1; icand < ic.size; ++icand) + { + irb = &(_m[ic.angle_index1[icand]][ic.angle_index2[icand]]); + wc = ic.weight[icand]; + //slog << "update from " << _m[ic.angle_index1[icand]][ic.angle_index2[icand]].fname << " " << wc << endl; + lc = irb->lbuf; + rc = irb->rbuf; + for (size_t i = 0; i < _ir_length; ++i) + { + left[i] += wc * lc[i]; + right[i] += wc * rc[i]; + } + } + // slog << "endof candidate update" << endl; +} + + +// convert an azimuth/elevation couple expressed in +// vertical polar coordinates to interaural polar coordinates +void HrtfCont::vertpol2interaurpol(float& az, float& el) const +{ + const float raz = az * DEG2RAD; + const float rel = el * DEG2RAD; + const float cosaz = cos(raz); + const float sinaz = sin(raz); + const float cosel = cos(rel); + const float sinel = sin(rel); + const float cosazcosel = cosaz * cosel; + + el = atan2(sinel, cosazcosel) * RAD2DEG; + az = atan2(cosel * sinaz, sqrt(cosazcosel * cosazcosel + sinel*sinel)) * RAD2DEG; +} diff --git a/src/hrtfcont.hpp b/src/hrtfcont.hpp new file mode 100644 index 0000000..304f99b --- /dev/null +++ b/src/hrtfcont.hpp @@ -0,0 +1,113 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#ifndef HRTFCONT_HPP_ +# define HRTFCONT_HPP_ + +# include <map> + +# include "ir_key.hpp" +# include "interpolation_candidates.hpp" + +// atomic element +// contains a buffer for a right and left impulse response +// corresponding to a given azimuth and elevation +typedef struct +{ + float *lbuf; // left buffer + float *rbuf; // right buffer + // string fname; +} ir_buffer; + + +// each HRTF is indexed by two angles in degree: azimuth and elevation +// Depending on the HRTF dabase considered, those angles may be expressed +// in a different coordinate system +// +// For the Listen Database, azimuth and elevation are expressed +// in the Vertical-polar coordinate system (most common spherical coordinate system) +// azimuth is in the range [0,360[ and distance between 2 azimuth is circular +// elevation is in the range [-90,+90] +// Consequently, database expressed in Vertical-polar coordinate system +// will be indexed by elevation first, and then azimuth +// +// For the CIPIC Database, azimuth and elevation are expressed +// in the Interaural-polar coordinate system +// elevation is in the range [0,360[ and distance between 2 elevations is circular +// azimuth is in the range [-90,+90] +// Consequently, database expressed in Interoral-polar coordinate system +// will be indexed by elevation first, and then azimuth + +// store impulse responses indexed by an angle whose values are in [0,360[ degree +// azimuth in vertical-polar coords, or elevation in interaural-polar coords +// corresponds to the second index of the storage map +typedef map<float, ir_buffer> angle2_cont; +// const interator on the second angle index +typedef angle2_cont::const_iterator angle2_cit; +// type used to store HRTFs corresponding of a given subject +// indexed by an angle in [-90,+90[ first and then by an angle in [0,360[ +typedef map<float, angle2_cont> hrtf_map; +// const interator on the first angle angle index +typedef map<float, angle2_cont>::const_iterator angle1_cit; + + +// HRTF Container +class HrtfCont +{ +public: + HrtfCont(const ir_key& k); + + // return the map indexing all available impulse responses + map<float, map<float, ir_buffer> >* map_get() {return &_m;} + + inline size_t ir_length_get() const {return _ir_length;} + + // update and interpolation structure containing the candidates to be interpolated + void set_candidates(interp_cdts& ic, float az, float el); + // called from set_candidates, for a fixed angle index 1 + // add the angle index2 candidates + void add_a2_candidates(interp_cdts& ic, float a1_key, float a2, float weight); + + // update filter coefficients from interpolation candidates + void update_from_candidates(const interp_cdts& ic, float* left, float* right); + +protected: + // normalize vertical polar coordinates expressed in degree + // ie: set elevation in [-90, +90], and azimuth in [-180, 180[ + void normalize_vertpolar_coords(float& azimuth, float& elevation) const; + + // angular distance (ie dist(350,0) == 10) + float angular_dist(float a1, float a2) const; + + // convert an azimuth/elevation couple expressed in + // vertical polar coordinates to interaural polar coordinates + void vertpol2interaurpol(float& az, float& el) const; + + size_t _ir_length; + hrtf_map _m; + const bool _vert_pol_coords; +}; + +#endif diff --git a/src/interpolation_candidates.cpp b/src/interpolation_candidates.cpp new file mode 100644 index 0000000..e54818b --- /dev/null +++ b/src/interpolation_candidates.cpp @@ -0,0 +1,37 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include "interpolation_candidates.hpp" + +// add a candidate represented by its indexes a1,a1 +// assumption: the method won't be called more than 4 times! +void interp_cdts::add(float a1, float a2, float w) +{ + angle_index1[size] = a1; + angle_index2[size] = a2; + weight[size] = w; + ++size; +} + diff --git a/src/interpolation_candidates.hpp b/src/interpolation_candidates.hpp new file mode 100644 index 0000000..c6c3334 --- /dev/null +++ b/src/interpolation_candidates.hpp @@ -0,0 +1,55 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#ifndef INTERPOLATION_CANDIDATES_HPP_ +# define INTERPOLATION_CANDIDATES_HPP_ + +# include <stddef.h> + +// this class is aimed to store the set of available measures, +// referenced by their azimuth and elevation +// required to interpolate the impulse response for a given position +// each candidate has an associated corresponding to its +// contribution in the final result +// a longer class name could be interpolation_candidates +class interp_cdts +{ +public: + interp_cdts() : size(0) {} + + void add(float a1, float a2, float w); + + // number of elements + size_t size; + // first angle index in the map in [-90,90] + float angle_index1[4]; + // second angle index in the map in [0,360[ + float angle_index2[4]; + // weights + float weight[4]; +}; + + +#endif diff --git a/src/ir_key.cpp b/src/ir_key.cpp new file mode 100644 index 0000000..caa3c96 --- /dev/null +++ b/src/ir_key.cpp @@ -0,0 +1,79 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include "ir_key.hpp" + +// The impulse response key comparison function +// returns true if its first argument is less than its second argument, +// and false otherwise +bool ir_key_isless::operator() (const ir_key& a, const ir_key& b) +{ + int tmp; + + if (a.minp_ap_dec < b.minp_ap_dec) + return true; + if (a.minp_ap_dec > b.minp_ap_dec) + return false; + + if (a.spectral < b.spectral) + return true; + if (a.spectral > b.spectral) + return false; + + if (a.length < b.length) + return true; + if (a.length > b.length) + return false; + + if (a.spectral < b.spectral) + return true; + if (a.spectral > b.spectral) + return false; + + tmp = a.path.compare(b.path); + if (tmp < 0) + return true; + if (tmp > 0) + return false; + + tmp = a.iir_regex.compare(b.iir_regex); + if (tmp < 0) + return true; + if (tmp > 0) + return false; + + if (a.is_azimuth_first < b.is_azimuth_first) + return true; + if (a.is_azimuth_first > b.is_azimuth_first) + return false; + + if (a.vertical_polar_coords < b.vertical_polar_coords) + return true; + if (a.vertical_polar_coords > b.vertical_polar_coords) + return false; + + // keys are equal, then !(a < b) + return false; +} diff --git a/src/ir_key.hpp b/src/ir_key.hpp new file mode 100644 index 0000000..45cffb3 --- /dev/null +++ b/src/ir_key.hpp @@ -0,0 +1,87 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#ifndef IR_KEY_HPP_ +# define IR_KEY_HPP_ + +# include <string> + +using namespace std; + +// this structure is aimed at identifying +// a requested impulse response set +// that will be use as it + +typedef struct s_ir_key +{ + // TODO + // * add ITD computation method!! + // * add resampling options + + // tells if we want a minphase/allpass decomposition + // having this options set to true implies more computations + // and a better hrtf interpolation + bool minp_ap_dec; + // tells whether we want to use the impulse response + // in the frequency domain or in the time domain + bool spectral; + // the length of the considered impulse response + size_t length; + // the path to the directory storing the set of impulse responses + string path; + // the regex that should match the name of valid IR files + // contains 2 groups corresponding to the azimuth and elevation + // composing the sound file name + string iir_regex; + // set to true if the azimuth is in the first group of iir_regex + // false if the elevation is first + bool is_azimuth_first; + // set to true in the db is indexed in vertical-polar coordinates + // set to false if the db is indexed in interaural-polar_coordinates + bool vertical_polar_coords; + + // constructor init the structure fields to default values + s_ir_key(): minp_ap_dec(false), + spectral(false), + length(0), + path(""), + iir_regex(""), + is_azimuth_first(true), + vertical_polar_coords(true) {} +} ir_key; + + +// The impulse response key comparison function +// returns true if its first argument is less than its second argument, +// and false otherwise +// it is necessary to index ir_key into sorted containers (map, etc...) +class ir_key_isless +{ +public: + bool operator() (const ir_key& a, const ir_key& b); +}; + + +#endif diff --git a/src/itdcont.cpp b/src/itdcont.cpp new file mode 100644 index 0000000..74a679b --- /dev/null +++ b/src/itdcont.cpp @@ -0,0 +1,92 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include "itdcont.hpp" +#include "flyweight_ir_factory.hpp" + +#include "logstring.hpp" + +// returns an itd set corresponding to key k +// assumption: constructor will alwas be called from FlyweightIrFactory +// and the field spectral set to false in the factory +ItdCont::ItdCont(const ir_key& k) +{ + ir_key k2 = k; + k2.minp_ap_dec = false; + + HrtfCont* hc = FlyweightIrFactory::instance()->hrtf_set_get(k2); + + for (angle1_cit ia1 = hc->map_get()->begin(); ia1 != hc->map_get()->end(); ia1++) + for (angle2_cit ia2 = ia1->second.begin(); ia2 != ia1->second.end(); ia2++) + _m[ia1->first][ia2->first] + = cross_correlation_itd(ia2->second.lbuf, ia2->second.rbuf, hc->ir_length_get()); +} + +// compute the cross-correlation between signals s1, s2 of size len_sig +// for n sample time lag +float ItdCont::cross_correlation(const float* s1, const float* s2, size_t len_sig, int n) +{ + float res = 0; + // slog << "ARG " << n << endl; + for (int i = (n >= 0 ? 0 : -n); i < (n <= 0 ? (int) len_sig : (int) len_sig - n); ++i) + res += s1[i] * s2[i+n]; + + // slog << "local cross corr " << n << ":" << res << endl; + return res; +} + +// compute an estimation of the ITD using crosscorrelation +// taking into account the left and right impulse responses +float ItdCont::cross_correlation_itd(const float *lsig, const float *rsig, size_t len_sig) +{ + int arg, argmax; + float val, valmax; + int ilen_sig = (int) len_sig; + + // slog << "lensig " << len_sig << ilen_sig<< endl; + + //slog << lsig << rsig << len_sig << endl; + argmax = 1-len_sig; + valmax = cross_correlation(lsig, rsig, len_sig, argmax); + for (arg = 2 - ilen_sig; arg < ilen_sig - 1; ++arg) + //{ + // slog << "FOR ARG " << arg << endl; + if ((val = cross_correlation(lsig, rsig, len_sig, arg)) > valmax) + { + argmax = arg; + valmax = val; + } + // } + // slog << "ARGMAXXXX " << argmax << endl; + return argmax; +} + +float ItdCont::itd_from_candidates(const interp_cdts& ic) +{ + float res = _m[ic.angle_index1[0] ][ic.angle_index2[0] ] * ic.weight[0]; + for (size_t i = 1; i < ic.size; ++i) + res += _m[ic.angle_index1[i]][ic.angle_index2[i]] * ic.weight[i]; + return res; +} diff --git a/src/itdcont.hpp b/src/itdcont.hpp new file mode 100644 index 0000000..1f6d36b --- /dev/null +++ b/src/itdcont.hpp @@ -0,0 +1,51 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#ifndef ITDCONT_HPP_ +# define ITDCONT_HPP_ + +# include <map> + +# include "ir_key.hpp" +# include "interpolation_candidates.hpp" + +class ItdCont +{ +public: + ItdCont(const ir_key& k); + // a positive return value means the right signal is late + // a negative return value means the left signal is late + float itd_from_candidates(const interp_cdts& ic); + + // to remove + inline map<float,map<float,float> >& map_get() {return _m;} +protected: + float cross_correlation_itd(const float *lsig, const float *rsig, size_t len_sig); + float cross_correlation(const float* s1, const float* s2, size_t len_sig, int n); + + map<float,map<float,float> > _m; +}; + +#endif diff --git a/src/logstring.cpp b/src/logstring.cpp new file mode 100644 index 0000000..8d4f5d0 --- /dev/null +++ b/src/logstring.cpp @@ -0,0 +1,43 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include "logstring.hpp" +#ifdef PD +# include <m_pd.h> +#endif + +stringstream slog; +const unsigned buffsize = 1024*32; +char cbuff[buffsize]; + +#ifdef PD +void logstring2pdterm() +{ + while(slog.getline(cbuff, buffsize)) + post(cbuff); + slog.clear(); +} + +#endif diff --git a/src/logstring.hpp b/src/logstring.hpp new file mode 100644 index 0000000..19e5447 --- /dev/null +++ b/src/logstring.hpp @@ -0,0 +1,42 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#ifndef LOGSTRING_HPP_ +# define LOGSTRING_HPP_ + +# include <iostream> +# include <sstream> + +using namespace std; + +// This string stream is aimed at storing cw_binaural~ instances messages +// in order to avoid the use of "post" which is PD dependant +extern stringstream slog; + +# ifdef PD +void logstring2pdterm(); +# endif + +#endif diff --git a/src/minphase_hrtfcont.cpp b/src/minphase_hrtfcont.cpp new file mode 100644 index 0000000..0b163dc --- /dev/null +++ b/src/minphase_hrtfcont.cpp @@ -0,0 +1,142 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include <math.h> +#include <m_pd.h> //needed for fft implementation +#include "minphase_hrtfcont.hpp" +#include "flyweight_ir_factory.hpp" + +MinPhaseHrtfCont::MinPhaseHrtfCont(const ir_key&k): + HrtfCont(k) +{ + ir_key k2 = k; + + // set impulse response length + _ir_length = k.length; + // FIXME: IR size verification ??? => power of 2 + + // get the container storing the coresponding temporal impulse response + k2.minp_ap_dec = false; + HrtfCont* hc = FlyweightIrFactory::instance()->hrtf_set_get(k2); + + // iterate on the container storing the temporal impulse response + for (angle1_cit ie = hc->map_get()->begin(); ie != hc->map_get()->end(); ie++) + for (angle2_cit ia2 = ie->second.begin(); ia2 != ie->second.end(); ia2++) + { + // current azimuth and elevation + const float el = ie->first; + const float az = ia2->first; + + // allocate buffers + float *lbuf = _m[el][az].lbuf = new float[_ir_length]; + float *rbuf = _m[el][az].rbuf = new float[_ir_length]; + + minphase_ir(ia2->second.lbuf, lbuf, _ir_length); + minphase_ir(ia2->second.rbuf, rbuf, _ir_length); + } +} + +// store in dst the magnitude of spectrum src of size n +// assumption: src is the specetrum of a real signal +void MinPhaseHrtfCont::magnitude(const float* src, float* dst, size_t n) +{ + for (size_t i = 1; i < n/2; ++i) + { + dst[i] = sqrt(src[i]*src[i] + src[n-i] * src[n-i]); + dst[n -i] = dst[i]; + } + dst[0] = src[0] >= 0 ? src[0] : -src[0]; + dst[n/2] = src[n/2] >= 0 ? src[n/2] : -src[n/2]; +} + +// store in dst the imaginary part of hilbert transform +// applied to signal src of size n +void MinPhaseHrtfCont::im_hilbert(const float* src, float* dst, size_t n) +{ + size_t i; + float tmp; + + if (dst != src) + for (i = 0; i < n; ++i) + dst[i] = src[i]; + + mayer_realfft(n, dst); + + for (i = 1; i < n/2; ++i) + { + tmp = dst[i]; + dst[i] = -dst[n -i]; + dst[n-i] = tmp; + } + dst[0] = dst[n/2] = 0; + + mayer_realifft(n, dst); + + for (i = 0; i < n; ++i) + dst[i] /= n; +} + +// store in dst the minphase impulse response corresponding to src +void MinPhaseHrtfCont::minphase_ir(const float* src, float* dst, int n) +{ + float *sig_spectrum = new float[n]; + float *magn, *phase; + int i; + + // compute the spectrum from input signal + for (i = 0; i < n; ++i) + sig_spectrum[i] = src[i]; + mayer_realfft(n, sig_spectrum); + + // get the magnitude of the spectrum => symetric signal + magn = new float[n]; + magnitude(sig_spectrum, magn, n); + + // compute hilbert transform of the log of the magnitude + // this computation correspond to the minimum phase + phase = new float[n]; + for (i = 0; i < n; ++i) + phase[i] = -log(magn[i]); + im_hilbert(phase, phase, n); + + float* real = new float[n]; + float* imag = new float[n]; + + for (i = 0; i < n; ++i) + { + real[i] = magn[i] * cos(phase[i]); + imag[i] = magn[i] * sin(phase[i]); + } + mayer_ifft(n, real, imag); + + delete [] sig_spectrum; + delete [] imag; + delete [] magn; + delete [] phase; + + for (i = 0; i < n; ++i) + dst[i] = real[i] /= n; + delete [] real; +} diff --git a/src/minphase_hrtfcont.hpp b/src/minphase_hrtfcont.hpp new file mode 100644 index 0000000..aed44e7 --- /dev/null +++ b/src/minphase_hrtfcont.hpp @@ -0,0 +1,45 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + + +#ifndef MINPHASE_HRTFCONT_HPP_ +# define MINPHASE_HRTFCONT_HPP_ + +# include "hrtfcont.hpp" + +// This class is used to store the minimal phase impulse response +// of a give HRTF +// These filters should be used in cunjunction with a pure delay +class MinPhaseHrtfCont: public HrtfCont +{ +public: + MinPhaseHrtfCont(const ir_key& k); +private: + void magnitude(const float* src, float* dst, size_t n); + void im_hilbert(const float* src, float* dst, size_t n); + void minphase_ir(const float* src, float* dst, int n); +}; + +#endif diff --git a/src/raw_wav_hrtfcont.cpp b/src/raw_wav_hrtfcont.cpp new file mode 100644 index 0000000..52bf82b --- /dev/null +++ b/src/raw_wav_hrtfcont.cpp @@ -0,0 +1,178 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include <regex.h> +#include <stdlib.h> +#include <dirent.h> +#include <sndfile.h> +#include <assert.h> +#include "raw_wav_hrtfcont.hpp" +#include "logstring.hpp" + + + +void RawWavHrtfCont::set_ir_buffer_from_path(ir_buffer& irb, string path) +{ + SF_INFO sfinfo; + SNDFILE* sndfd; + + // slog << "Trying " << path << std::endl; + //irb.fname = path.c_str(); + + // opening the sound file + sfinfo.format = 0; + sndfd = sf_open(path.c_str(), SFM_READ, &sfinfo); + if (!sndfd) + { + slog << "failed to open " << path << endl; + throw 0; + } + + //slog << sfinfo.channels << endl; + if (sfinfo.channels != 2) + { + slog << path << " has " << sfinfo.channels << " channels instead of 2!" << endl; + throw 0; + } + assert(sfinfo.samplerate == 44100); + + //slog << "frames" << sfinfo.frames << endl; + const size_t wav_length = sfinfo.frames;// / 2; + //slog << "ir length before" << _ir_length << endl; + if (!_ir_length) + // we're extrating the first impulse response + // all other impulse response extrated + // must have the same length + _ir_length = wav_length; + //slog << "ir length" << _ir_length << endl; + + // if the considered file is empty + // or if its size is different from previously + // extracted impulse response files + assert(wav_length && (wav_length == _ir_length)); + + float* wavdata = new float[sfinfo.frames*2]; + sf_read_float(sndfd, wavdata, sfinfo.frames*2); + + irb.lbuf = new float[wav_length]; + irb.rbuf = new float[wav_length]; + + + for (size_t i = 0; i < wav_length; ++i) + { + irb.lbuf[i] = wavdata[i*2]; + irb.rbuf[i] = wavdata[i*2+1]; + } + + delete [] wavdata; + sf_close(sndfd); +} + +// constructor should parse the directory +// deduce the how many different azimuths +// and elevations are present, according +// to a given regexp +// => intermediate list structure +// then use a faster structure to store +// the extracted wavs!!! +// FIXME: check the map size is != 0 +RawWavHrtfCont::RawWavHrtfCont(const ir_key& k): + HrtfCont(k) +{ + // regex vars + int re_status; + regex_t re; + regmatch_t pmatch[3]; + // directory parsing vars + DIR *dp; + struct dirent *dirp; + // other + float az, elev; + + _ir_length = 0; + + // slog << "Raw Wav HRTF CONT " << k.iir_regex << endl; + // compile the regex + if (regcomp(&re, k.iir_regex.c_str(), REG_EXTENDED|REG_ICASE)) + { + slog << "regex: " << k.iir_regex << " is not supported!" << endl; + throw 0; + } + + // test existence of directory supposed to store hrtf impulse responses + if((dp = opendir(k.path.c_str())) == NULL) + { + slog << "directory " << k.path << " does not exist!" << endl; + throw 0; + } + + // extract wav files matching the regex + // contained in path + while ((dirp = readdir(dp))) + //{ slog << dirp->d_name << endl; + if (!regexec(&re, dirp->d_name, 3, pmatch, REG_NOTBOL|REG_NOTEOL)) + { + //slog << "match!" << endl; + string dirname(dirp->d_name); + string group0 = dirname.substr(pmatch[1].rm_so, pmatch[1].rm_eo); + string group1 = dirname.substr(pmatch[2].rm_so, pmatch[2].rm_eo); + + if (k.is_azimuth_first) + { + az = atof(group0.c_str()); + elev = atof(group1.c_str()); + } + else + { + az = atof(group1.c_str()); + elev = atof(group0.c_str()); + } + // slog << "(az,elev)=(" << az << ", " << elev << ") for file " << dirname << endl; + + if (k.vertical_polar_coords) + { + normalize_vertpolar_coords(az, elev); + set_ir_buffer_from_path(_m[elev][az], k.path + "/" + dirname); + } + else + { + // TO BE IMPROVED + if (elev > 180) + elev -= 360; + set_ir_buffer_from_path(_m[az][elev], k.path + "/" + dirname); + } + + } + + // close the directory structure + closedir(dp); + + // free the regex + regfree(&re); + + //slog << "END OF RAW WAV HRTF CONT" << endl; +} + + diff --git a/src/raw_wav_hrtfcont.hpp b/src/raw_wav_hrtfcont.hpp new file mode 100644 index 0000000..8748202 --- /dev/null +++ b/src/raw_wav_hrtfcont.hpp @@ -0,0 +1,53 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#ifndef RAW_WAV_HRTFCONT_HPP_ +# define RAW_WAV_HRTFCONT_HPP_ + +# include "hrtfcont.hpp" + +// to remove + +# include <string> + +class RawWavHrtfCont: public HrtfCont +{ + // constructor should parse the directory + // deduce the how many different azimuths + // and elevations are present, according + // to a given regexp + // => intermediate list structure + // then use a faster structure to store + // the extracted wavs!!! +public: + RawWavHrtfCont(const ir_key& k); + +private: + void set_ir_buffer_from_path(ir_buffer& irb, string path); + +}; + + +#endif diff --git a/src/riff.cpp b/src/riff.cpp new file mode 100644 index 0000000..80bc87b --- /dev/null +++ b/src/riff.cpp @@ -0,0 +1,60 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + +#include "riff.hpp" + + +Riff::Riff(int size): GenericRiff(size) +{ + _fifo = new float[size]; + _fifo_id = size-1; + for (int i = 0; i < size; ++i) + _fifo[i] = 0; +} + +Riff::~Riff() +{ + + delete [] _fifo; +} + +void Riff::process(float* input, float* output, int n) +{ + int i,j; + float* curcoeff; + + for (i = 0; i < n; ++i) + { + float out = 0; + _fifo[_fifo_id] = input[i]; + + for (j = _fifo_id, curcoeff = _coeffs; j < _size; ++j, ++curcoeff) + out += _fifo[j] * *curcoeff; + for (j = 0; j < _fifo_id; ++j, ++curcoeff) + out += _fifo[j] * *curcoeff; + + output[i] = out; + _fifo_id = _fifo_id ? _fifo_id - 1 : _size - 1; + } +} diff --git a/src/riff.hpp b/src/riff.hpp new file mode 100644 index 0000000..0d34245 --- /dev/null +++ b/src/riff.hpp @@ -0,0 +1,43 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + + +#ifndef RIFF_HPP_ +# define RIFF_HPP_ + +# include "generic_riff.hpp" + +class Riff: public GenericRiff +{ +public: + Riff(int size); + ~Riff(); + virtual void process(float* input, float* output, int n); +protected: + int _fifo_id; + float *_fifo; +}; + +#endif diff --git a/src/spectral_hrtfcont.cpp b/src/spectral_hrtfcont.cpp new file mode 100644 index 0000000..0710506 --- /dev/null +++ b/src/spectral_hrtfcont.cpp @@ -0,0 +1,67 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include <m_pd.h> //needed for fft implementation +#include "spectral_hrtfcont.hpp" +#include "flyweight_ir_factory.hpp" + +SpectralHrtfCont::SpectralHrtfCont(const ir_key& k): + HrtfCont(k) +{ + ir_key k2 = k; + + // length is twice the impulse response considered + _ir_length = k.length * 2; + // FIXME: IR size verification ??? => power of 2 + + // get the container storing the coresponding temporal impulse response + k2.spectral = false; + HrtfCont* hc = FlyweightIrFactory::instance()->hrtf_set_get(k2); + + // iterate on the container storing the temporal impulse response + for (angle1_cit ie = hc->map_get()->begin(); ie != hc->map_get()->end(); ie++) + for (angle2_cit ia2 = ie->second.begin(); ia2 != ie->second.end(); ia2++) + { + // current azimuth and elevation + const float el = ie->first; + const float az = ia2->first; + + // allocate buffers + float *lbuf = _m[el][az].lbuf = new float[_ir_length]; + float *rbuf = _m[el][az].rbuf = new float[_ir_length]; + + // fill buffers + for (size_t i = 0; i < k.length; ++i) + { + lbuf[i] = ia2->second.lbuf[i]; + rbuf[i] = ia2->second.rbuf[i]; + lbuf[i+k.length] = rbuf[i+k.length] = 0; + } + + // apply fft to the buffers + mayer_realfft(_ir_length, lbuf); + mayer_realfft(_ir_length, rbuf); + } +} diff --git a/src/spectral_hrtfcont.hpp b/src/spectral_hrtfcont.hpp new file mode 100644 index 0000000..422fb4f --- /dev/null +++ b/src/spectral_hrtfcont.hpp @@ -0,0 +1,43 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + + +#ifndef SPECTRAL_HRTFCONT_HPP_ +# define SPECTRAL_HRTFCONT_HPP_ + +# include "hrtfcont.hpp" + +// This class is used to store spectral domain impulse responses +// The corresponding temporal impulse response of size n +// will be 0-padded to size 2n before putting it to spectral domain +// In order to be able to use the spectral impulse responses +// in an overlap-add filter +class SpectralHrtfCont: public HrtfCont +{ +public: + SpectralHrtfCont(const ir_key& k); +}; + +#endif diff --git a/src/unprocessed_fixed_size_hrtfcont.cpp b/src/unprocessed_fixed_size_hrtfcont.cpp new file mode 100644 index 0000000..08691b0 --- /dev/null +++ b/src/unprocessed_fixed_size_hrtfcont.cpp @@ -0,0 +1,70 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + +#include "unprocessed_fixed_size_hrtfcont.hpp" +#include "flyweight_ir_factory.hpp" + +UnprocessedFixedSizeHrtfCont::UnprocessedFixedSizeHrtfCont(const ir_key& k): + HrtfCont(k) +{ + HrtfCont* hc; + ir_key k2 = k; + + // set the container's length + _ir_length = k.length; + + // get the container containing the raw extracted wav + k2.length = 0; + hc = FlyweightIrFactory::instance()->hrtf_set_get(k2); + + // iterate on the container storing the raw wav impulse response + for (angle1_cit ie = hc->map_get()->begin(); ie != hc->map_get()->end(); ie++) + for (angle2_cit ia2 = ie->second.begin(); ia2 != ie->second.end(); ia2++) + { + // get current azimuth and elevation + const float el = ie->first; + const float az = ia2->first; + + // allocate buffers + _m[el][az].lbuf = new float[k.length]; + _m[el][az].rbuf = new float[k.length]; + + // copy available data + const size_t av_data_len = _ir_length < hc->ir_length_get() ? _ir_length : hc->ir_length_get(); + for (size_t i = 0; i < av_data_len; ++i) + { + _m[el][az].lbuf[i] = ia2->second.lbuf[i]; + _m[el][az].rbuf[i] = ia2->second.rbuf[i]; + } + + // padd with 0 missing data + for (size_t i = av_data_len; i < _ir_length; ++i) + { + _m[el][az].lbuf[i] = 0; + _m[el][az].rbuf[i] = 0; + } + + } +} diff --git a/src/unprocessed_fixed_size_hrtfcont.hpp b/src/unprocessed_fixed_size_hrtfcont.hpp new file mode 100644 index 0000000..ec127c8 --- /dev/null +++ b/src/unprocessed_fixed_size_hrtfcont.hpp @@ -0,0 +1,42 @@ +/* + cw_binaural~: a binaural synthesis external for pure data + by David Doukhan - david.doukhan@gmail.com - http://www.limsi.fr/Individu/doukhan + and Anne Sedes - sedes.anne@gmail.com + Copyright (C) 2009-2011 David Doukhan and Anne Sedes + + For more details, see CW_binaural~, a binaural synthesis external for Pure Data + David Doukhan and Anne Sedes, PDCON09 + + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. +*/ + + + +#ifndef UNPROCESSED_FIXED_SIZE_HRTFCONT_HPP_ +# define UNPROCESSED_FIXED_SIZE_HRTFCONT_HPP_ + +# include "hrtfcont.hpp" + +// This container store an unprocessed impulse response +// of size fixed in the constructory argument +// if the requested size is bigger than the available impulse +// response it is padded with 0 +class UnprocessedFixedSizeHrtfCont: public HrtfCont +{ +public: + UnprocessedFixedSizeHrtfCont(const ir_key& k); +}; + +#endif |
