diff general/numerical/frint.m @ 4:e44f49929e56

Adding reorganised general toolbox, now in several subdirectories.

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/general/numerical/frint.m	Sat Jan 12 19:21:22 2013 +0000
@@ -0,0 +1,22 @@
+function y=frint(r,x,dc)
+% frint - Fractional integration
+%
+% y=frint(r,x[,dc])
+% Returns the rth order integral of x, where r can be fractional
+% if dc=0 (default is 1) dc term is not added back in
+
+x=x(:);
+n=length(x);
+z=fft(x); z0=z(1); z(1)=0;
+K=(0:n-1)';
+D=exp(2*pi*i*K/n)-1;
+D(1)=1; % to avoid divide by zero
+z2=z./D.^r;
+y=real(ifft(z2)); % integral minus DC term
+
+% put DC term  back as a ramp
+if nargin<3 || dc==1,
+	y=y+sign(z0)*(K*abs(z0)/n).^r;
+end
+
+