#include #include #include /********************************************************************/ /* fonction de fit "f(x;h,x0,sig,c) = h*exp(-0.5*((x-x0)/sig)^2)+c" */ /* avec nvar=1 npar=4 */ /* x = x[0] */ /* p[0]=h, p[1]=x0, p[2]=sig, p[3]=c */ /********************************************************************/ double gauss1(double const* x,double const* p) { double xc; xc = (x[0]-p[1])/p[2]; xc *= -xc/2.; return((xc>=-80.) ? p[0]*exp(xc)+p[3] : p[3]); } /* la definition de cette fonction n'est pas obligatoire */ /* dp[0] = df(x;h,x0,sig,c)/dh, dp[1] = etc... */ double gauss1_der(double const* x,double const* p,double* dp) { double xc,xc2,e,f; xc = (x[0]-p[1])/p[2]; xc2 = xc*xc; e = -xc2/2.; if(e>=-80.) e = exp(e); else e = 0.; f = p[0]*e; dp[0] = e; dp[1] = xc / p[2] *f; dp[2] = xc2 / p[2] *f; dp[3] = 1.; return(f+p[3]); } /********************************************************/ /* fonction de fit gaussienne 2D avec correlation */ /* f(x,y;h,x0,y0,sx,sy,rho,c) = */ /* h * exp{ -0.5*[ (x-x0)^2/sx^2 + (y-y0)^2/s^2 */ /* -2*rho*(x-x0)*(y-y0)/(sx*sy) ]} */ /* avec nvar=2 npar=7 */ /* x = x[0] et y=x[2] */ /* p[0]=h, p[1]=x0, p[2]=y0 */ /* p[3]=sx, p[4]=sy, p[5]=rho, p[6]=fond */ /********************************************************/ double gauss2(double const* x,double const* p) { double ax,ay,a; ax = (x[0]-p[1])/p[3]; ay = (x[1]-p[2])/p[4]; a = -(ax*ax+ay*ay-2.*p[5]*ax*ay)/2.; return((a>=-80.) ? p[0]*exp(a)+p[6] : p[6]); } /* la definition de cette fonction n'est pas obligatoire */ /* dp[0] = f(x,y;h,x0,y0,sx,sy,rho,c)/dh, dp[1] = etc... */ double gauss2_der(double const* x,double const* p,double* dp) { double ax,ay,a,p0ea; ax = (x[0]-p[1])/p[3]; ay = (x[1]-p[2])/p[4]; a = -(ax*ax+ay*ay-2.*p[5]*ax*ay)/2.; p0ea = (a>=-80.) ? exp(a) : 0.; dp[0] = p0ea; p0ea *= p[0]; dp[1] = (ax-p[5]*ay)/p[3] *p0ea; dp[2] = (ay-p[5]*ax)/p[4] *p0ea; dp[3] = ax*dp[1]; dp[4] = ay*dp[2]; dp[5] = ax*ay *p0ea; dp[6] = 1.; return(p0ea+p[6]); }