(*Input values(all masses in MeV)*)
q[M_] =1/2*Sqrt[(M^2-4*Mpi^2)];

a[M_,G_] =M^2-12*Mpi^2+(8*Pi*q[M]^3)/G+(24*q[M]*Mpi^2/M)*Log[(M+2*q[M])/(2*Mpi)];
b[M_,G_] =a[M,G]*(M^2+2*Mpi^2)/(6*Mpi^2*M^2)-2*q[M]^2*(G*M^2+8*Pi*q[M]^3)/(3*G*Mpi^2*M^2);

Mpi=139.5;

Mr=775.1;
Gr=150.4;
mr=0.819*10^6;
ar=a[Mr,Gr];
br=b[Mr,Gr];

Mr1=1406;
Gr1=455;
mr1=0.217*10^6;
ar1=a[Mr1,Gr1];
br1=b[Mr1,Gr1];

Mr2=1700;
Gr2=235;
mr2=150*10^6;
ar2=a[Mr2,Gr2];
br2=b[Mr2,Gr2];

Na=1.03;
Be=-0.121;
Ga=0.032;

(*range of s*)
smin=4*Mpi^2;
smax=4*10^6;

(*defining functions*)
p[s_] =(s-4*Mpi^2);
sigma[s_] =Sqrt[(p[s]/s)];
tau[s_] =Sqrt[(-p[s]/s)];

j[s_] =p[s]*(sigma[s]*Log[(sigma[s]-1)/(sigma[s]+1)]+2);
fbar[s_] =(1+s/mr^2)*(a[Mr,Gr]/(a[Mr,Gr]-b[Mr,Gr]*s-j[s]));
fbar1[s_] =(1+s/mr1^2)*(a[Mr1,Gr1]/(a[Mr1,Gr1]-b[Mr1,Gr1]*s-j[s]));
fbar2[s_] =(1+s/mr2^2)*(a[Mr2,Gr2]/(a[Mr2,Gr2]-b[Mr2,Gr2]*s-j[s]));
(*Main Function*)
Fgs[s_] =Na*(fbar[s]+Be*fbar1[s]+Ga*fbar2[s])/(1+Be+Ga);

Plot[Norm[Fgs[EE^2]],{EE,280,1300}]

Plot[Arg[Fgs[EE^2]](*+2 Pi HeavisideTheta[EE-1220]-2 Pi HeavisideTheta[EE-1582]*),{EE,0,2000}]

Plot[Arg[Fgs[EE^2]]+2 Pi HeavisideTheta[EE-1220]-2 Pi HeavisideTheta[EE-1582],{EE,0,2000}]

Limit[Fgs[s],s->0]

Plot[Im[Fgs[EE^2]],{EE,0,1300}]

Plot[Re[Fgs[EE^2]],{EE,0,1300}]

Clear[smax]

smax=200000^2;

Omnes[s_]:=Exp[s/Pi (NIntegrate[(Arg[Fgs[x]]+2 Pi HeavisideTheta[x-1220^2]- 2 Pi HeavisideTheta[x-1582^2])/(x(x-s)),{x, 4 Mpi^2, s-10^(-7)}]+NIntegrate[(Arg[Fgs[x]]+2 Pi HeavisideTheta[x-1220^2]-2 Pi HeavisideTheta[x-1582^2])/(x(x-s)),{x,s+10^(-7),smax}]-0*2.77/Pi*(( π Log[1-(s/smax)])/s))]

Plot[Im[Fgs[EE^2]],{EE,1570,1590}]

Omnes[1000^2]

Table[Abs[Fgs[(200+n 100)^2]],{n,0,10}]

Table[Abs[Omnes[(200+n 100)^2]],{n,0,10}]

Table[Abs[Fgs[(200+n 100)^2]],{n,0,10}]/1.03

delta00[s_]:=Arg[Fgs[s]]180/Pi

k2[s_]=Sqrt[s-4 mK^2]/2; mK=0.495*1000;Mpr=1.5*1000; e1=2.4; e2=-5.5;
eta00inp[s_]=1-(e1 k2[s]/Sqrt[s]+e2 k2[s]^2/s)(Mpr^2-s)/s;

eta00[s_]:=1 HeavisideTheta[4 mK^2-s]+eta00inp[s] HeavisideTheta[s-4 mK^2]

mpi=139;

Ref0[s_]:=Re[(eta00[s] Exp[2 I delta00[s] Pi/180]-1)/(2 I Sqrt[(s-4 mpi^2)/s])]

Imf0[s_]:=Im[(eta00[s] Exp[2 I delta00[s] Pi/180]-1)/(2 I Sqrt[(s-4 mpi^2)/s])]

f0[s_]:=(Ref0[s]+I Imf0[s])Sqrt[(s-4 mpi^2)/s]

ParametricPlot[{Re@#,Im@#}&@f0[EE^2],{EE,300,1400},PlotRange->All,PlotStyle->Thick,Epilog->{Red,PointSize[0.03]}]

Plot[eta00[EE^2],{EE,300,1400}]

Plot[delta00[EE^2],{EE,300,1400}]