annotate aim-mat/modules/usermodule/mellin/CalMI_Rich.old.and.slow.m @ 4:537f939baef0 tip

various bug fixes and changed copyright message
author Stefan Bleeck <bleeck@gmail.com>
date Tue, 16 Aug 2011 14:37:17 +0100
parents 20ada0af3d7d
children
rev   line source
tomwalters@0 1 % tester generating function for Calculation of the Mellin Image in 'aim-mat'
tomwalters@0 2 %
tomwalters@0 3 % INPUT VALUES: SAI3d : 3D SAI
tomwalters@0 4 % NAPparam: Parameter for NAP
tomwalters@0 5 % SAIparam: Parameter for SAI
tomwalters@0 6 % MIparam: Parameter for MI
tomwalters@0 7 % RETURN VALUE:
tomwalters@0 8 % MI3d: 3D Mellin Image
tomwalters@0 9 %
tomwalters@0 10 % (c) 2003, University of Cambridge, Medical Research Council
tomwalters@0 11 % Original Code IRINO T
tomwalters@0 12 % 10 Jan. 2002
tomwalters@0 13 % Modified by R. Turner (ret26@cam.ac.uk)
tomwalters@0 14 % Feb. 2003
bleeck@3 15 % (c) 2011, University of Southampton
bleeck@3 16 % Maintained by Stefan Bleeck (bleeck@gmail.com)
bleeck@3 17 % download of current version is on the soundsoftware site:
bleeck@3 18 % http://code.soundsoftware.ac.uk/projects/aimmat
bleeck@3 19 % documentation and everything is on http://www.acousticscale.org
bleeck@3 20
tomwalters@0 21
tomwalters@0 22 function [MI3d] = CalMI_Rich(SAI3d,options,sample_rate)
tomwalters@0 23
tomwalters@0 24 %here we setup all the variables that the function requires
tomwalters@0 25 %%% CHANGE THE NAMES SO IT'S CLEAR WHAT'S GOING ON!!!
tomwalters@0 26
tomwalters@0 27 %what do these do?
tomwalters@0 28 MIparam.NSAIPhsCmp = 0; %was 2 but aim uses a window starting at 0ms
tomwalters@0 29 MIparam.F0mode = 300;
tomwalters@0 30
tomwalters@0 31 %these values set the resolution and scale of the axes in the final MT
tomwalters@0 32 %they are now defined in the parameter file
tomwalters@0 33 MIparam.TFval = options.TFval;
tomwalters@0 34 MIparam.c_2pi = options.c_2pi;
tomwalters@0 35 Lenc2pi = length(MIparam.c_2pi);
tomwalters@0 36 LenTF = length(MIparam.TFval);
tomwalters@0 37
tomwalters@0 38 %not sure what these do
tomwalters@0 39 MIparam.Mu = -0.5; % flat if Mu <0.5: high pass, Mu>0.5 low pass
tomwalters@0 40 SAIparam.Nwidth = 0; %sets the negative width of the window
tomwalters@0 41 MIparam.SSI = options.ssi;
tomwalters@0 42
tomwalters@0 43 %we find the maximum temporal interval size, the no of channels and the no
tomwalters@0 44 %of frames
tomwalters@0 45 [NumCh, LenSAI, LenFrame] = size(SAI3d);
tomwalters@0 46
tomwalters@0 47 %going to remove any information in the first N channels
tomwalters@0 48 % N=5;
tomwalters@0 49 % for no_ch = 1:N,
tomwalters@0 50 % SAI3d(no_ch,:,:)=zeros(LenSAI,LenFrame);
tomwalters@0 51 % end;
tomwalters@0 52
tomwalters@0 53 %saving the SAI3d for comparison with irino's
tomwalters@0 54 % savefile = 'SAI3d_AIM.mat';
tomwalters@0 55 % SAI3d_AIM=SAI3d/max(max(max(SAI3d)));
tomwalters@0 56 % save(savefile,'SAI3d_AIM');
tomwalters@0 57
tomwalters@0 58 %for using irino's SAI in AIM
tomwalters@0 59 % load SAI3d_Irino.mat;
tomwalters@0 60 % SAI3d=SAI3d_Irino;
tomwalters@0 61
tomwalters@0 62 %here we set up the Frs values for each channel
tomwalters@0 63 cf_afb = [100 6000]; %4500]; altered for 2dAT
tomwalters@0 64 NAPparam.Frs = FcNch2EqERB(min(cf_afb),max(cf_afb),NumCh);
tomwalters@0 65 NAPparam.fs = sample_rate;
tomwalters@0 66
tomwalters@0 67 %We initialise the mellin image matrix (why is it best to do this - is it
tomwalters@0 68 %so we can spot errors more easily?
tomwalters@0 69 MI3d = zeros(Lenc2pi,LenTF,LenFrame);
tomwalters@0 70
tomwalters@0 71 %user information (commented out) + new wait bar
tomwalters@0 72 %disp('*** MI Calculation ***');
tomwalters@0 73
tomwalters@0 74 %MIparam.RangeAudFig = [];
tomwalters@0 75 %disp(MIparam);
tomwalters@0 76 waithand=waitbar(0,'generating the MT');
tomwalters@0 77
tomwalters@0 78 % set the frame range, the mellin image is calculated for
tomwalters@0 79 if (options.do_all_frames == 1)
tomwalters@0 80 start_frame = 1;
tomwalters@0 81 end_frame = LenFrame;
tomwalters@0 82 else
tomwalters@0 83 start_frame = options.framerange(1);
tomwalters@0 84 end_frame = options.framerange(2);
tomwalters@0 85 end;
tomwalters@0 86
tomwalters@0 87
tomwalters@0 88 %this section does the filter response alignment
tomwalters@0 89 for nfr = start_frame:end_frame
tomwalters@0 90 %set up the waitbar
tomwalters@0 91 fraction_complete=nfr/LenFrame;
tomwalters@0 92 waitbar(fraction_complete);
tomwalters@0 93
tomwalters@0 94 %generate the new matricies of the frames
tomwalters@0 95 SAIval = SAI3d(:,:,nfr);
tomwalters@0 96 SAIPhsCmp = zeros(size(SAIval));
tomwalters@0 97
tomwalters@0 98 %we shift each channel along the time interval axis by adding zeros
tomwalters@0 99 for nch = 1:NumCh,
tomwalters@0 100 NPeriod = NAPparam.fs/NAPparam.Frs(nch) * MIparam.NSAIPhsCmp;
tomwalters@0 101 shift_matrix = [zeros(1,fix(NPeriod)), SAIval(nch,:)];
tomwalters@0 102 SAIPhsCmp(nch,1:LenSAI) = shift_matrix(1:LenSAI);
tomwalters@0 103 %MI3d(:,:,nfr) = SAIPhsCmp'; %added line
tomwalters@0 104 %SAIPhsCmp = SAI3d(:,:,nfr); %tester line
tomwalters@0 105 end;
tomwalters@0 106 if MIparam.F0mode == 0
tomwalters@0 107 % No estimation of F0
tomwalters@0 108 else
tomwalters@0 109 F0est(nfr) = MIparam.F0mode;
tomwalters@0 110 end;
tomwalters@0 111
tomwalters@0 112 %Here we extract the information required to ensure that we have
tomwalters@0 113 %only one presentation of the 'timbre' information
tomwalters@0 114 ZeroLoc = abs(SAIparam.Nwidth)*NAPparam.fs/1000+1;
tomwalters@0 115 MarginAF = 0; % No margin by introducing WinAF
tomwalters@0 116
tomwalters@0 117 % maah: set the range for the auditory image
tomwalters@0 118 if (options.do_all_image == 1)
tomwalters@0 119 MIparam.RangeAudFig = [1 LenSAI];
tomwalters@0 120 else
tomwalters@0 121 MIparam.RangeAudFig = options.audiorange;
tomwalters@0 122 end;
tomwalters@0 123 % maah: was MIparam.RangeAudFig = [ZeroLoc+[0 , (fix(NAPparam.fs/F0est(nfr))-MarginAF)]];
tomwalters@0 124
tomwalters@0 125 % Calculation of the Mellin Coefficients
tomwalters@0 126 [MImtrx] = CalMellinCoef_Rich(SAIPhsCmp,NAPparam,MIparam);
tomwalters@0 127
tomwalters@0 128 %Output into the 3d matrix
tomwalters@0 129 MI3d(1:Lenc2pi,1:LenTF,nfr) = MImtrx;
tomwalters@0 130
tomwalters@0 131 % maah: Magnitude
tomwalters@0 132 MI3d(1:Lenc2pi,1:LenTF,nfr) = abs(MI3d(1:Lenc2pi,1:LenTF,nfr));
tomwalters@0 133 %MI3d(:,:,nfr) = MImtrx; %added line
tomwalters@0 134
tomwalters@0 135 end;
tomwalters@0 136
tomwalters@0 137 close(waithand);