Skip to content

Instantly share code, notes, and snippets.

@DanHickstein
Created September 17, 2018 01:02
Show Gist options
  • Save DanHickstein/230be11b8dff66d913e5888fc978be31 to your computer and use it in GitHub Desktop.
Save DanHickstein/230be11b8dff66d913e5888fc978be31 to your computer and use it in GitHub Desktop.
clear all;
N=1001;
Rm = fix(N/2)+1;
NBF=500;
I=(1:N)';
R2 = ((I-Rm).^2);
R = I-Rm;
M = zeros(N,NBF);
Mc =zeros(N,NBF);
Mc(:,1)=exp(-R2);
M(:,1)=2*exp(-R2);
for n=1:(NBF-1)
angn=exp(n^2-2*n^2*log(n)+sum(log(0.5:(n^2-0.5))));
M(Rm,n+1)=2*angn;
for r=1:(N-Rm)
Mc(r+Rm,n+1)=exp(n^2-r^2+2*n^2*log(r/n));
Mc(Rm-r,n+1)=Mc(r+Rm,n+1);
aux=2*Mc(r+Rm,n+1)+2*angn*Mc(r+Rm,1);
for l=max(1,r^2-100):min(n^2-1,r^2+100)
aux=aux+2*exp(n^2-r^2-n^2*log(n^2)+l*log(r^2)+sum(log((l+1):n^2))+sum(log(0.5:(n^2-l-0.5)))-sum(log(1:(n^2-l))));
end;
M(r+Rm,n+1)=aux;
M(Rm-r,n+1)=aux;
end;
n
end;
save('e:\basis1000_1.bst','Mc', '-ASCII');
save('e:\basis1000pr_1.bst','M', '-ASCII');
return
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment