File:MagaplotFragment.jpg

Fragment of the image http://mizugadro.mydns.jp/t/index.php/File:Magaplot300.jpg

Explicit plot of function maga and related functions.

Thin green curve: $~ y=\,$ nori$(x)$ $\,=\,$ mori$\big(\sqrt{x}\big)^2\,=\,$ $\displaystyle \frac{J_0\big(L\, \sqrt{x}\big)^2}{(1\!-\!x)^2}$

where $J_0=\,$BesselJ0, the zeroth Bessel function, and $L$ is its first zero, $L=\,$BesselJZero$[0,1]$.

The thick blue and thick red lines represent the half of the Bessel transform of square of mori with quadratic phase;

naga$(x) = \displaystyle 2 \,\int_0^\infty \, $mori$(p )^2 \, \exp(\mathrm i \,x \,p^2) \, p \mathrm d p$ $\, =\displaystyle \int_0^\infty \, $nori$(q) \, \exp(\mathrm i \,x \,q)\, \mathrm d q$

The thick black curve shows the maga function,

$y=\,$maga$(x) \, = 1-|$naga$(x)|^2$

The dotted curve shows its fit,

$y=\mathrm{fit2}(x)=\displaystyle \frac{ x^{3/2}}{\sqrt{5-x/2+x^2/4+x^3}}$

The thin blue line shows deviation of this fit from the maga function; in order to mage this deviation visible, it is scaled with factor 50:

$y=50\big( \mathrm{fit2}(x)-\mathrm{maga}(x)\big)$

C++ generator of curves
//using namespace std;
 * 1) include
 * 2) include
 * 3) include 
 * 4) include
 * 1) define z_type std::complex
 * 2) define Re(x) x.real
 * 3) define Im(x) x.imag
 * 4) define RI(x) x.real,x.imag
 * 5) define DB double
 * 6) define DO(x,y) for(x=0;x<y;x++)

void fft(z_type *a, int N, int o) // Arry is FIRST! {z_type u,w,t; int i,j,k,l,e=1,L,p,q,m; q=N/2; p=2; for(m=1;p<N;m++) p*=2; p=N-1; z_type y=z_type(0.,M_PI*o); j=0; for(i=0;i<p;i++){ if(i<j) { t=a[j]; a[j]=a[i]; a[i]=t;} k=q; while(k<=j){ j-=k; k/=2;} j+=k; } for(l=0;l<m;l++) { L=e; e*=2; u=1.; w=exp(y/((double)L)); for(j=0;j<L;j++) { for(i=j;i<N;i+=e){k=i+L; t=a[k]*u; a[k]=a[i]-t; a[i]+=t;} u*=w; } } }

void ado(FILE *O, int X, int Y) {      fprintf(O,"%c!PS-Adobe-2.0 EPSF-2.0\n",'%'); fprintf(O,"%c%cBoundingBox: 0 0 %d %d\n",'%','%',X,Y); fprintf(O,"/M {moveto} bind def\n"); fprintf(O,"/L {lineto} bind def\n"); fprintf(O,"/S {stroke} bind def\n"); fprintf(O,"/s {show newpath} bind def\n"); fprintf(O,"/C {closepath} bind def\n"); fprintf(O,"/F {fill} bind def\n"); fprintf(O,"/o {.01 0 360 arc C S} bind def\n"); fprintf(O,"/times-Roman findfont 20 scalefont setfont\n"); fprintf(O,"/W {setlinewidth} bind def\n"); fprintf(O,"/RGB {setrgbcolor} bind def\n");}

// #include "ado.cin" // #include"fafo.cin" //#include "mori.cin" DB L1= 2.404825557695773; DB L2= 5.5200781102863115; DB L3= 8.653727912911013; DB morin(DB x){ return j0(L1*x)/(1-x*x);} DB mori0(DB x){ int n,m; DB s, xx=x*x; DB c[16]={ 1., -0.4457964907366961303, 0.07678538241994023453, -0.0071642885058902232688, 0.00042159522055140947688, -0.000017110542281627483109, 5.0832583976057607495e-7, -1.1537378620148452816e-8, 2.0662789231930073316e-10, -2.9948657413756059965e-12, 3.5852738451127332173e-14,-3.6050239634659700777e-16, 3.0877184831292878827e-18, -2.2798156440952688462e-20, 1.4660907878585489441e-22,-8.2852774398657968065e-25}; //  16th term seems to fail; perhaps, due to the C++ rounding errors. //with m=15, at |x|<2, the error is of order of 10^(-16) //In this sense, the result is accurate while |x|<2. m=15; s=c[m]*xx; for(n=m-1;n>0;n--){ s+=c[n]; s*=xx;} return 1.+s;} DB mori(DB x){if(fabs(x)<2) return mori0(x); return morin(x);}

DB fit(DB x){ DB a,b,c,cc,s;    a=-.055; b=0.02; c=0.45; cc=c*c; s=x*(a+x*(b+x*cc)); return c*x*sqrt(x/(1.+s)); }

DB fit2(DB x){ DB q; q=5.+x*(-.5+x*(.25+x)); return x*sqrt(x/q);}

int main{z_type * a, *b, c; int j,m,n,N; FILE *o; DB scale,step,x,y,q,p,u,v; n=15; N=pow(2,n); //scale=.5; step=sqrt(2*M_PI/N); scale=100*step; DB dx=step*scale; DB dp=step/scale; printf("2^%2d=%8d scale=%6.4lf dx=%9.8lf dp=%9.8lf\n",n,N,scale,dx,dp); a=(z_type *) malloc((size_t)((N+1)*sizeof(z_type))); b=(z_type *) malloc((size_t)((N+1)*sizeof(z_type)));

DO(m,N){x=m*dx; y=mori(sqrt(x)); y*=y; a[m]=y; b[m]=y*dx; } b[0]*=.5; fft(b,N,1); //DO(j,N)  printf("%2d %18.15f %18.15f  %18.15f %18.15f\n", j, RI(a[j]), RI(b[j])  ); //o=fopen("32.eps","w"); ado(o,1008,228); o=fopen("magaplo.eps","w"); ado(o,1008,208); fprintf(o,"4 104 translate 100 100 scale\n 1 setlinejoin 2 setlinecap\n"); for(m=0;m<11;m++) {M(m,-1) L(m,1)} for(n=-1;n<2;n++) {M(0,n) L(10,n)} fprintf(o,".002 W S\n");
 * 1) define M(x,y) fprintf(o,"%6.4f %6.4f M\n",0.+x,0.+y);
 * 2) define L(x,y) fprintf(o,"%6.4f %6.4f L\n",0.+x,0.+y);
 * 3) define o(x,y) fprintf(o,"%6.4f %6.4f o\n",0.+x,0.+y);

DO(m,N){x=dx*m; y=Re(a[m]); if(m==0)M(x,y)else L(x,y);if(x>10) break;} fprintf(o,".008 W 0 .8 0 RGB S\n"); //DO(n,N){p=dp*n; y=fit(p); if(n==0)M(p,y)else L(p,y);if(p>10) break;} fprintf(o,".008 W 1 0 0 RGB S\n"); //DO(n,N){p=dp*n; y=fit2(p); if(n==0)M(p,y)else L(p,y);if(p>10) break;} fprintf(o,".008 W 0 0 1 RGB S\n"); fprintf(o,".01 W 0 0 1 RGB\n"); for(n=4;n10) break;}

DO(n,N){p=dp*n; y=Re(b[n]); if(n==0)M(p,y)else L(p,y);if(p>10) break;} fprintf(o,".02 W 0 0 1 RGB S\n"); DO(n,N){p=dp*n; y=Im(b[n]); if(n==0)M(p,y)else L(p,y);if(p>10) break;} fprintf(o,".02 W 1 0 0 RGB S\n"); DO(n,N){p=dp*n; u=Re(b[n]); v=Im(b[n]); y=1-u*u-v*v; if(n==0)M(p,y)else L(p,y);if(p>10) break;} fprintf(o,".02 W 0 0 0 RGB S\n"); // fprintf(o,".01 W 0 0 0 RGB\n"); //for(n=4;n10) break;}

//DO(n,N){p=dp*n; u=Re(b[n]); v=Im(b[n]); y=fit(p)-(1.-u*u-v*v); y*=50; if(n==0)M(p,y)else L(p,y);if(p>10) break;} fprintf(o,".004 W .8 0 0 RGB S\n"); DO(n,N){p=dp*n; u=Re(b[n]); v=Im(b[n]); y=fit2(p)-(1.-u*u-v*v); y*=50; if(n==0)M(p,y)else L(p,y);if(p>10) break;} fprintf(o,".008 W 0 0 1 RGB S\n");

fprintf(o,"showpage\n%cTrailer",'%'); fclose(o); system("epstopdf magaplo.eps"); system(   "open magaplo.pdf"); free(a); free(b); }

Latex generator of labels
\documentclass[12pt]{article} \usepackage{geometry} \usepackage{graphics} \usepackage{rotating} \paperwidth 420pt \paperheight 216pt \textwidth 1420pt \textheight 300pt \topmargin -108pt \oddsidemargin -73pt \newcommand \ds {\displaystyle} \newcommand \sx {\scalebox} \newcommand \rme {\mathrm{e}} \newcommand \rot {\begin{rotate}} \newcommand \ero {\end{rotate}} \newcommand \ing {\includegraphics} \pagestyle{empty} \parindent 0pt \begin{document} \begin{picture}(410,214) \put(10,0){\ing{magaplo}} %\put(.6,120){\sx{1.25}{$y$}} \put(2,199){\sx{1.3}{$1$}} %\put(.6,56){\sx{1.4}{$\frac{1}{2}$}} \put(2, 100){\sx{1.3}{$0$}} \put(-2, 0){\sx{1.3}{$-\!1$}} \put(111,88){\sx{1.33}{$1$}} \put(211,88){\sx{1.33}{$2$}} \put(311,89){\sx{1.33}{$3$}} \put(406,89){\sx{1.4}{$x$}} %\put(411,89){\sx{1.33}{$4$}} %\put(512,89){\sx{1.33}{$5$}} %\put(613,89){\sx{1.33}{$6$}} %\put(713,89){\sx{1.33}{$7$}} %\put(813,89){\sx{1.33}{$8$}} %\put(914,89){\sx{1.33}{$9$}} %\put(1010,89){\sx{1.4}{$x$}} %\put(15,192){\rot{-36}\sx{1.3}{$y\!=\! \mathrm{nori}(x)$}\ero}

\put(280,199){\rot{3}\sx{1.3}{$y\!=\! \mathrm{fit2}(x))$}\ero} \put(280,178){\rot{3}\sx{1.3}{$y\!=\! \mathrm{maga}(x))$}\ero}

\put(118,130){\rot{-16}\sx{1.3}{$y\!=\! \mathrm{nori}(x)$}\ero} \put(292,142){\rot{-4}\sx{1.3}{$y\!=\! \Im(\mathrm{naga}(x))$}\ero} \put(292,120){\rot{-4}\sx{1.3}{$y\!=\! \Re(\mathrm{naga}(x))$}\ero}

\put(170,26){\sx{1.3}{$y\!=\! 50 \Big(\mathrm{fit2}(x)-\mathrm{maga}(x)\Big)$}} \end{picture}\end{document}