Index: trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h
===================================================================
--- trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h	(revision 2834)
+++ trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h	(revision 2835)
@@ -132,6 +132,6 @@
   BlindPixelFitFunc fgSinglePheFitFunc;     //! In the beginning, 
   Int_t             fgSinglePheFitNPar;     //! we want to be flexible using different functions
-  
-  inline static Double_t gfKto4(Double_t *x, Double_t *par)
+
+  inline static Double_t fAna(Double_t *x, Double_t *par)
     {
 
@@ -152,5 +152,4 @@
       if (sigma1 < sigma0)
         return NoWay;
-      
       
       Double_t mu2 = (2.*mu1)-mu0;  
@@ -188,7 +187,64 @@
       return TMath::Exp(-1.*lambda)*par[5]*sum;
       
-    } //->
-
-  inline static Double_t gfKto5(Double_t *x, Double_t *par)
+    }
+    
+  inline static Double_t fKto4(Double_t *x, Double_t *par)
+    {
+
+      Double_t lambda = par[0];  
+      
+      Double_t sum = 0.;
+      Double_t arg = 0.;
+      
+      Double_t mu0 = par[1];
+      Double_t mu1 = par[2];
+      
+      if (mu1 < mu0)
+        return NoWay;
+
+      Double_t sigma0 = par[3];
+      Double_t sigma1 = par[4];
+      
+      if (sigma1 < sigma0)
+        return NoWay;
+      
+      Double_t mu2 = (2.*mu1)-mu0;  
+      Double_t mu3 = (3.*mu1)-(2.*mu0);
+      Double_t mu4 = (4.*mu1)-(3.*mu0);
+      
+      Double_t sigma2 = TMath::Sqrt((2.*sigma1*sigma1) - (sigma0*sigma0));  
+      Double_t sigma3 = TMath::Sqrt((3.*sigma1*sigma1) - (2.*sigma0*sigma0));
+      Double_t sigma4 = TMath::Sqrt((4.*sigma1*sigma1) - (3.*sigma0*sigma0));
+      
+      Double_t lambda2 = lambda*lambda;
+      Double_t lambda3 = lambda2*lambda;
+      Double_t lambda4 = lambda3*lambda;
+      
+      // k=0:
+      arg = (x[0] - mu0)/sigma0;
+      sum = TMath::Exp(-0.5*arg*arg)/sigma0;
+      
+      // k=1:
+      arg = (x[0] - mu1)/sigma1;
+      sum += lambda*TMath::Exp(-0.5*arg*arg)/sigma1;
+      
+      // k=2:
+      arg = (x[0] - mu2)/sigma2;
+      sum += 0.5*lambda2*TMath::Exp(-0.5*arg*arg)/sigma2;
+      
+      // k=3:
+      arg = (x[0] - mu3)/sigma3;
+      sum += 0.1666666667*lambda3*TMath::Exp(-0.5*arg*arg)/sigma3;
+      
+      // k=4:
+      arg = (x[0] - mu4)/sigma4;
+      sum += 0.041666666666667*lambda4*TMath::Exp(-0.5*arg*arg)/sigma4;
+      
+      return TMath::Exp(-1.*lambda)*par[5]*sum;
+      
+    } 
+
+  
+  inline static Double_t fKto5(Double_t *x, Double_t *par)
     {
       
@@ -255,5 +311,5 @@
   
   
-  inline static Double_t gfKto6(Double_t *x, Double_t *par)
+  inline static Double_t fKto6(Double_t *x, Double_t *par)
     {
       
