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