annotate aim-mat/tools/@signal/generatedampsinus.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
bleeck@3 1 % method of class @signal
bleeck@3 2 % function sig=generatedampsinus(sig,carfre,modfre,amplitude,halflife)
bleeck@3 3 % INPUT VALUES:
bleeck@3 4 % sig: original @signal with length and samplerate
bleeck@3 5 % carfre: carrier frequency (Hz) [1000]
bleeck@3 6 % modfre: modulation frequency (Hz) [100]
bleeck@3 7 % amplitude: [1]
bleeck@3 8 % halflife: time for the envelope envelope to decrease exponentielly
bleeck@3 9 % to 1/2
bleeck@3 10 %
bleeck@3 11 %
bleeck@3 12 % RETURN VALUE:
bleeck@3 13 % sig: @signal
bleeck@3 14 %
bleeck@3 15 % This external file is included as part of the 'aim-mat' distribution package
bleeck@3 16 % (c) 2011, University of Southampton
bleeck@3 17 % Maintained by Stefan Bleeck (bleeck@gmail.com)
bleeck@3 18 % download of current version is on the soundsoftware site:
bleeck@3 19 % http://code.soundsoftware.ac.uk/projects/aimmat
bleeck@3 20 % documentation and everything is on http://www.acousticscale.org
bleeck@3 21
bleeck@3 22
bleeck@3 23
bleeck@3 24 function sig=generatedampsinus(sig,carfre,modfre,amplitude,halflife,jitter)
bleeck@3 25 % generates a damped sinusoid, that is a carrier pure tone modulated with a
bleeck@3 26 % exponentially decreasing envelope.
bleeck@3 27 % sig is the signal
bleeck@3 28 % carfre is the carrier frequency of the pure tone
bleeck@3 29 % modfre is the modulation frequency (in Hz)
bleeck@3 30 % amplitude is the final amplitude
bleeck@3 31 % halflife is the time in seconds, in which the envelope drops to its half value
bleeck@3 32 % if halflife is 'gamma', instead a gammaenvelope is used
bleeck@3 33
bleeck@3 34 % jitter is how regular the pulses are (0-1) 1=100%
bleeck@3 35 % if jitter is 'repeat', only one pulse is generated and repeated
bleeck@3 36
bleeck@3 37
bleeck@3 38 if nargin < 6
bleeck@3 39 jitter=0;
bleeck@3 40 end
bleeck@3 41 if nargin < 5
bleeck@3 42 halflife=0.01;
bleeck@3 43 end
bleeck@3 44 if nargin < 4
bleeck@3 45 amplitude=1;
bleeck@3 46 end
bleeck@3 47
bleeck@3 48 if nargin < 3
bleeck@3 49 modfre=100;
bleeck@3 50 end
bleeck@3 51 if nargin < 2
bleeck@3 52 carfre=1000;
bleeck@3 53 end
bleeck@3 54
bleeck@3 55
bleeck@3 56
bleeck@3 57
bleeck@3 58 if isequal(halflife,'gamma') && isequal(jitter,'repeat')
bleeck@3 59 % code for Tom's 2010 MSc project
bleeck@3 60
bleeck@3 61 % generate only one period and repeat it n times
bleeck@3 62 sr=getsr(sig);
bleeck@3 63 periodlen=1/modfre;
bleeck@3 64 oneperiod=signal(periodlen,sr);
bleeck@3 65
bleeck@3 66 sinus=generatesinus(oneperiod,carfre,amplitude,0);
bleeck@3 67
bleeck@3 68 % calculate envelope and mult both
bleeck@3 69 envelope=oneperiod;
bleeck@3 70 t=[0:1/sr:periodlen-1/sr]*500*(500/getnrpoints(envelope));
bleeck@3 71 env=power(t,4).*exp(-2*pi*t/2);
bleeck@3 72 env=env./max(env);
bleeck@3 73 envelope=setvalues(envelope,env);
bleeck@3 74 envelope=envelope/max(envelope)*amplitude;
bleeck@3 75 envelope=setstarttime(envelope,0);
bleeck@3 76
bleeck@3 77 % set the envelope and the amplitude
bleeck@3 78 oneperiod=sinus*envelope;
bleeck@3 79
bleeck@3 80 % repeat for the right number of repeats
bleeck@3 81 nr_repeats=round(getlength(sig)/getlength(envelope));
bleeck@3 82 fsig=oneperiod;
bleeck@3 83
bleeck@3 84 if nr_repeats>1
bleeck@3 85 for i=1:nr_repeats-1
bleeck@3 86 fsig=append(fsig,oneperiod);
bleeck@3 87 end
bleeck@3 88 end
bleeck@3 89
bleeck@3 90 sig=fsig;
bleeck@3 91
bleeck@3 92 sig=setname(sig,sprintf('damped sinusoid %4.2f kHz, Modulation=%4.1f Hz, gamma envelope',carfre/1000,modfre));
bleeck@3 93
bleeck@3 94 % figure(1)
bleeck@3 95 % plot(env)
bleeck@3 96
bleeck@3 97 return
bleeck@3 98 end
bleeck@3 99
bleeck@3 100
bleeck@3 101
bleeck@3 102
bleeck@3 103
bleeck@3 104
bleeck@3 105
bleeck@3 106
bleeck@3 107
bleeck@3 108
bleeck@3 109
bleeck@3 110
bleeck@3 111
bleeck@3 112
bleeck@3 113
bleeck@3 114
bleeck@3 115
bleeck@3 116
bleeck@3 117 sinus=generatesinus(sig,carfre,amplitude,0);
bleeck@3 118
bleeck@3 119 % calculate envelope and mult both
bleeck@3 120 envelope=sig;
bleeck@3 121
bleeck@3 122 time_const=halflife/0.69314718;
bleeck@3 123
bleeck@3 124 env_vals=getvalues(envelope);
bleeck@3 125 sr=getsr(envelope);
bleeck@3 126 reprate=1/modfre;
bleeck@3 127
bleeck@3 128
bleeck@3 129 sig_len=getlength(sig);
bleeck@3 130 % when regular, all modulataions are at the same time:
bleeck@3 131 if jitter==0
bleeck@3 132 pulse_times=0:reprate:sig_len;
bleeck@3 133 else
bleeck@3 134 nr_pulses=ceil(sig_len/reprate);
bleeck@3 135 modulation_period=1/modfre;
bleeck@3 136 pulse_times(1)=0;
bleeck@3 137 for n = 2:nr_pulses
bleeck@3 138 jittering=(rand-0.5*jitter)*modulation_period;
bleeck@3 139 pulse_times(n) = pulse_times(n-1) + modulation_period+jittering;
bleeck@3 140 end
bleeck@3 141 end
bleeck@3 142
bleeck@3 143
bleeck@3 144
bleeck@3 145 % oldval=0;
bleeck@3 146 % corr=exp(1/sr/time_const);
bleeck@3 147 % t=0;
bleeck@3 148 next_pulse=pulse_times(1);
bleeck@3 149 pulse_counter=1;
bleeck@3 150
bleeck@3 151 % for i=1:getnrpoints(envelope);
bleeck@3 152 % oldval=oldval/corr;
bleeck@3 153 % t=t+1/sr;
bleeck@3 154 % if t>next_pulse
bleeck@3 155 % oldval=1;
bleeck@3 156 % pulse_counter=pulse_counter+1;
bleeck@3 157 % if length(pulse_times)>=pulse_counter
bleeck@3 158 % next_pulse=pulse_times(pulse_counter);
bleeck@3 159 % else
bleeck@3 160 % next_pulse=inf;
bleeck@3 161 % end
bleeck@3 162 % end
bleeck@3 163 % env_vals(i)= oldval;
bleeck@3 164 % end
bleeck@3 165
bleeck@3 166 % onsettimeconstant
bleeck@3 167 %
bleeck@3 168 t=[0:1/sr:reprate-1/sr]*500;
bleeck@3 169 env=power(t,4).*exp(-2*pi*t/2);
bleeck@3 170 env=env./max(env);
bleeck@3 171
bleeck@3 172 % e=linspace(env(length(env)-l),0,l+1);
bleeck@3 173 % env(length(env)-l:end)=e;
bleeck@3 174
bleeck@3 175 % e=linspace(1,0,l+1);
bleeck@3 176 % env(length(env)-l:end)=env(length(env)-l:end).*e.*e;
bleeck@3 177
bleeck@3 178 t=0;
bleeck@3 179 ct=1;
bleeck@3 180 for i=1:getnrpoints(envelope);
bleeck@3 181 t=t+1/sr;
bleeck@3 182 ct=ct+1;
bleeck@3 183 if t>next_pulse
bleeck@3 184 ct=1;
bleeck@3 185 pulse_counter=pulse_counter+1;
bleeck@3 186 if length(pulse_times)>=pulse_counter
bleeck@3 187 next_pulse=pulse_times(pulse_counter);
bleeck@3 188 else
bleeck@3 189 next_pulse=inf;
bleeck@3 190 end
bleeck@3 191 end
bleeck@3 192 if ct>length(env)
bleeck@3 193 env_vals(i)= 0;
bleeck@3 194 else
bleeck@3 195 env_vals(i)= env(ct);
bleeck@3 196 end
bleeck@3 197 end
bleeck@3 198
bleeck@3 199
bleeck@3 200
bleeck@3 201 %
bleeck@3 202 % onsettime=0;
bleeck@3 203 % for i=1:getnrpoints(envelope);
bleeck@3 204 % t=t+1/sr;
bleeck@3 205 %
bleeck@3 206 % % env_vals(i)= exp(-(time)/time_const);
bleeck@3 207 % env_vals(i)= power(t,onsettime)*exp(-(t)/time_const);
bleeck@3 208 % t=mod(t,reprate);
bleeck@3 209 %
bleeck@3 210 % end
bleeck@3 211 envelope=setvalues(envelope,env_vals);
bleeck@3 212 envelope=envelope/max(envelope)*amplitude;
bleeck@3 213 envelope=setstarttime(envelope,0);
bleeck@3 214
bleeck@3 215 % set the envelope and the amplitude
bleeck@3 216 sig=sinus*envelope;
bleeck@3 217
bleeck@3 218 sig=setname(sig,sprintf('damped sinusoid %4.2f kHz, Modulation=%4.1f Hz, halflife=%4.1f ms',carfre/1000,modfre,halflife*1000));
bleeck@3 219
bleeck@3 220 % plot(sig)
bleeck@3 221
bleeck@3 222
bleeck@3 223
bleeck@3 224
bleeck@3 225
bleeck@3 226
bleeck@3 227
bleeck@3 228
bleeck@3 229
bleeck@3 230
bleeck@3 231
bleeck@3 232
bleeck@3 233
bleeck@3 234
bleeck@3 235
bleeck@3 236
bleeck@3 237 %
bleeck@3 238 %
bleeck@3 239 %
bleeck@3 240 % % from here on old code, might not work!
bleeck@3 241 %
bleeck@3 242 % sinus=generatesinus(sig,carfre,amplitude,0);
bleeck@3 243 %
bleeck@3 244 % % calculate envelope and mult both
bleeck@3 245 % envelope=sig;
bleeck@3 246 %
bleeck@3 247 % time_const=halflife/0.69314718;
bleeck@3 248 %
bleeck@3 249 % env_vals=getvalues(envelope);
bleeck@3 250 % sr=getsr(envelope);
bleeck@3 251 % reprate=1/modfre;
bleeck@3 252 %
bleeck@3 253 %
bleeck@3 254 % sig_len=getlength(sig);
bleeck@3 255 % % when regular, all modulataions are at the same time:
bleeck@3 256 % if jitter==0
bleeck@3 257 % pulse_times=0:reprate:sig_len;
bleeck@3 258 % else
bleeck@3 259 % nr_pulses=ceil(sig_len/reprate);
bleeck@3 260 % modulation_period=1/modfre;
bleeck@3 261 % pulse_times(1)=0;
bleeck@3 262 % for n = 2:nr_pulses
bleeck@3 263 % jittering=(rand-0.5*jitter)*modulation_period;
bleeck@3 264 % pulse_times(n) = pulse_times(n-1) + modulation_period+jittering;
bleeck@3 265 % end
bleeck@3 266 % end
bleeck@3 267 %
bleeck@3 268 %
bleeck@3 269 %
bleeck@3 270 % % oldval=0;
bleeck@3 271 % % corr=exp(1/sr/time_const);
bleeck@3 272 % % t=0;
bleeck@3 273 % next_pulse=pulse_times(1);
bleeck@3 274 % pulse_counter=1;
bleeck@3 275 %
bleeck@3 276 % % for i=1:getnrpoints(envelope);
bleeck@3 277 % % oldval=oldval/corr;
bleeck@3 278 % % t=t+1/sr;
bleeck@3 279 % % if t>next_pulse
bleeck@3 280 % % oldval=1;
bleeck@3 281 % % pulse_counter=pulse_counter+1;
bleeck@3 282 % % if length(pulse_times)>=pulse_counter
bleeck@3 283 % % next_pulse=pulse_times(pulse_counter);
bleeck@3 284 % % else
bleeck@3 285 % % next_pulse=inf;
bleeck@3 286 % % end
bleeck@3 287 % % end
bleeck@3 288 % % env_vals(i)= oldval;
bleeck@3 289 % % end
bleeck@3 290 %
bleeck@3 291 % % onsettimeconstant
bleeck@3 292 % %
bleeck@3 293 % t=[0:1/sr:reprate-1/sr]*500;
bleeck@3 294 % env=power(t,4).*exp(-2*pi*t/2);
bleeck@3 295 % env=env./max(env);
bleeck@3 296 %
bleeck@3 297 % % e=linspace(env(length(env)-l),0,l+1);
bleeck@3 298 % % env(length(env)-l:end)=e;
bleeck@3 299 %
bleeck@3 300 % % e=linspace(1,0,l+1);
bleeck@3 301 % % env(length(env)-l:end)=env(length(env)-l:end).*e.*e;
bleeck@3 302 %
bleeck@3 303 % t=0;
bleeck@3 304 % ct=1;
bleeck@3 305 % for i=1:getnrpoints(envelope);
bleeck@3 306 % t=t+1/sr;
bleeck@3 307 % ct=ct+1;
bleeck@3 308 % if t>next_pulse
bleeck@3 309 % ct=1;
bleeck@3 310 % pulse_counter=pulse_counter+1;
bleeck@3 311 % if length(pulse_times)>=pulse_counter
bleeck@3 312 % next_pulse=pulse_times(pulse_counter);
bleeck@3 313 % else
bleeck@3 314 % next_pulse=inf;
bleeck@3 315 % end
bleeck@3 316 % end
bleeck@3 317 % if ct>length(env)
bleeck@3 318 % env_vals(i)= 0;
bleeck@3 319 % else
bleeck@3 320 % env_vals(i)= env(ct);
bleeck@3 321 % end
bleeck@3 322 % end
bleeck@3 323 %
bleeck@3 324 %
bleeck@3 325 %
bleeck@3 326 % %
bleeck@3 327 % % onsettime=0;
bleeck@3 328 % % for i=1:getnrpoints(envelope);
bleeck@3 329 % % t=t+1/sr;
bleeck@3 330 % %
bleeck@3 331 % % % env_vals(i)= exp(-(time)/time_const);
bleeck@3 332 % % env_vals(i)= power(t,onsettime)*exp(-(t)/time_const);
bleeck@3 333 % % t=mod(t,reprate);
bleeck@3 334 % %
bleeck@3 335 % % end
bleeck@3 336 % envelope=setvalues(envelope,env_vals);
bleeck@3 337 % envelope=envelope/max(envelope)*amplitude;
bleeck@3 338 % envelope=setstarttime(envelope,0);
bleeck@3 339 %
bleeck@3 340 % % set the envelope and the amplitude
bleeck@3 341 % sig=sinus*envelope;
bleeck@3 342 %
bleeck@3 343 % sig=setname(sig,sprintf('damped sinusoid %4.2f kHz, Modulation=%4.1f Hz, halflife=%4.1f ms',carfre/1000,modfre,halflife*1000));
bleeck@3 344 %
bleeck@3 345 % % plot(sig)