comparison userProgramsTim/mellin_trafo.m @ 38:c2204b18f4a2 tip

End nov big change

+function out = mellin_trafo(inx,iny)
+
+%This function computes the Mellin Transformation of a one-dimensional
+%Signal in analytical terms given as
+%
+% S(p)=Integral from 0 to infty of s(t) t^(p-1)dt
+% Which equals 
+% S(c) = Integral from -Infty to Infty s(t)*exp(j*c*log(t))d(log(t))
+%
+% taken from Irino and Patterson, Speech Communication 2002 Eq. (1) and (3)
+% Tim Juergens, September 2011
+%
+
+%Resample the Signal onto a log(t) axis
+minimalx=min(inx);
+maximalx=max(inx);
+logarithmicx= exp([log(minimalx):(log(maximalx)-log(minimalx))/length(inx):log(maximalx)]);
+logarithmicy= interp1(inx,iny,logarithmicx,'linear','extrap');
+%figure, semilogx(logarithmicx,logarithmicy);
+
+%Absolute of the inverse Fourier-Transform of the resampled signal
+out = abs(ifft(logarithmicy));
+%figure, plot(out(1:40));
+
+
+%%just to show that the mellin transform results in invariable patterns if
+%%the formants are a constant ratio
+% frequencies_short = [1:8:8000];
+% frequencies_long = [1:10:10000];
+% Intervals_long = 1./frequencies_long;
+% Intervals_short = 1./frequencies_short;
+% contoursin = sin(2*pi*0.00015.*frequencies_long).^2
+% figure, plot(frequencies_long,contoursin), hold on, plot(frequencies_short,contoursin,'r')
+% xlabel('Frequency (Hz)')
+% ylabel('Rate (arb. units)')
+% figure, semilogx(Intervals_long,contoursin), hold on, plot(Intervals_short,contoursin,'r')
+% xlabel('Interval (s)')
+% ylabel('Rate (arb. units)')
+% m1 = mellin_trafo(Intervals_long,contoursin);
+% m2 = mellin_trafo(Intervals_short,contoursin);
+% figure, plot(m1(1:40)), hold on, plot(m2(1:40),'r');
+% xlabel('Mellin variable c');
+% ylabel('Magnitude');
+
+