]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ANALYSIS/AliUnfolding.cxx
Fixed grep syntax on OS X in the Analysis Plugin
[u/mrichter/AliRoot.git] / ANALYSIS / AliUnfolding.cxx
index 82d0acc09248bf17782b640813e3fdf6b2e6a4e3..ded255871637371a0865decae9d6db4c3647dd2a 100644 (file)
@@ -55,7 +55,7 @@ Float_t AliUnfolding::fgRegularizationWeight = 10000;
 Int_t AliUnfolding::fgSkipBinsBegin = 0;
 Float_t AliUnfolding::fgMinuitStepSize = 0.1;                 // (usually not needed to be changed) step size in minimization
 Float_t AliUnfolding::fgMinuitPrecision = 1e-6;               // minuit precision
-Int_t   AliUnfolding::fgMinuitMaxIterations = 1000000;            // minuit maximum number of iterations\r
+Int_t   AliUnfolding::fgMinuitMaxIterations = 1000000;            // minuit maximum number of iterations
 Double_t AliUnfolding::fgMinuitStrategy = 1.;                 // minuit strategy
 Bool_t AliUnfolding::fgMinimumInitialValue = kFALSE;          // set all initial values at least to the smallest value among the initial values
 Float_t AliUnfolding::fgMinimumInitialValueFix = -1;
@@ -197,7 +197,7 @@ void AliUnfolding::SetBayesianParameters(Float_t smoothing, Int_t nIterations)
   fgBayesianSmoothing = smoothing;
   fgBayesianIterations = nIterations;
 
-  Printf("AliUnfolding::SetBayesianParameters --> Parameters set to %d iterations with smoothing %f", fgBayesianIterations, fgBayesianSmoothing);\r
+  Printf("AliUnfolding::SetBayesianParameters --> Parameters set to %d iterations with smoothing %f", fgBayesianIterations, fgBayesianSmoothing);
 }
 
 //____________________________________________________________________
@@ -276,12 +276,20 @@ void AliUnfolding::SetStaticVariables(TH2* correlation, TH1* measured, TH1* effi
     fgEntropyAPriori = new TVectorD(fgMaxParams);
   if (!fgEfficiency)
     fgEfficiency = new TVectorD(fgMaxParams);
-  if (!fgUnfoldedAxis)
+  if (fgUnfoldedAxis)
+  {
     delete fgUnfoldedAxis;
-  fgUnfoldedAxis = new TAxis(*(correlation->GetXaxis()));
-  if (!fgMeasuredAxis)
+    fgUnfoldedAxis = 0;
+  }
+  if (!fgUnfoldedAxis) 
+    fgUnfoldedAxis = new TAxis(*(correlation->GetXaxis()));
+  if (fgMeasuredAxis)
+  {
     delete fgMeasuredAxis;
-  fgMeasuredAxis = new TAxis(*(correlation->GetYaxis()));    
+    fgMeasuredAxis = 0;
+  }
+  if (!fgMeasuredAxis) 
+    fgMeasuredAxis = new TAxis(*(correlation->GetYaxis()));
 
   fgCorrelationMatrix->Zero();
   fgCorrelationCovarianceMatrix->Zero();
@@ -295,14 +303,14 @@ void AliUnfolding::SetStaticVariables(TH2* correlation, TH1* measured, TH1* effi
     if (sum <= 0)
       continue;
     Float_t maxValue = 0;
-    Int_t maxBin = -1;
+    //    Int_t maxBin = -1;
     for (Int_t j=1; j<=correlation->GetNbinsY(); ++j)
     {
       // find most probably value
       if (maxValue < correlation->GetBinContent(i, j))
       {
         maxValue = correlation->GetBinContent(i, j);
-        maxBin = j;
+       //        maxBin = j;
       }
 
       // npart sum to 1
@@ -555,7 +563,7 @@ Int_t AliUnfolding::UnfoldWithBayesian(TH2* correlation, TH1* aEfficiency, TH1*
   Double_t* result = new Double_t[kMaxT];
   Double_t* efficiency = new Double_t[kMaxT];
   Double_t* binWidths = new Double_t[kMaxT];
-  Double_t* normResponse = new Double_t[kMaxT]; \r
+  Double_t* normResponse = new Double_t[kMaxT]; 
 
   Double_t** response = new Double_t*[kMaxT];
   Double_t** inverseResponse = new Double_t*[kMaxT];
@@ -563,7 +571,7 @@ Int_t AliUnfolding::UnfoldWithBayesian(TH2* correlation, TH1* aEfficiency, TH1*
   {
     response[i] = new Double_t[kMaxM];
     inverseResponse[i] = new Double_t[kMaxM];
-    normResponse[i] = 0;\r
+    normResponse[i] = 0;
   }
 
   // for normalization
@@ -577,7 +585,7 @@ Int_t AliUnfolding::UnfoldWithBayesian(TH2* correlation, TH1* aEfficiency, TH1*
     {
       response[t][m] = correlation->GetBinContent(t+1, m+1);
       inverseResponse[t][m] = 0;
-      normResponse[t] += correlation->GetBinContent(t+1, m+1);\r
+      normResponse[t] += correlation->GetBinContent(t+1, m+1);
     }
   }
 
@@ -593,13 +601,13 @@ Int_t AliUnfolding::UnfoldWithBayesian(TH2* correlation, TH1* aEfficiency, TH1*
     prior[t] = measuredCopy[t];
     result[t] = 0;
     binWidths[t] = aResult->GetXaxis()->GetBinWidth(t+1);
-\r
-    for (Int_t m=0; m<kMaxM; m++) { // Normalise response matrix\r
-      if (normResponse[t] != 0) \r
-       response[t][m] /= normResponse[t];\r
-      else\r
-        Printf("AliUnfolding::UnfoldWithBayesian: Empty row,column in response matrix, for truth bin %d",t);\r
-    }\r
+
+    for (Int_t m=0; m<kMaxM; m++) { // Normalise response matrix
+      if (normResponse[t] != 0) 
+       response[t][m] /= normResponse[t];
+      else
+        Printf("AliUnfolding::UnfoldWithBayesian: Empty row,column in response matrix, for truth bin %d",t);
+    }
   }
 
   // pick prior distribution
@@ -629,7 +637,7 @@ Int_t AliUnfolding::UnfoldWithBayesian(TH2* correlation, TH1* aEfficiency, TH1*
     {
       Float_t norm = 0;
       for (Int_t t = kStartBin; t<kMaxT; t++)
-        norm += response[t][m] * prior[t] * efficiency[t];\r
+        norm += response[t][m] * prior[t] * efficiency[t];
 
       // calc. chi2: (measured - response * prior) / error
       if (measuredError[m] > 0)
@@ -641,7 +649,7 @@ Int_t AliUnfolding::UnfoldWithBayesian(TH2* correlation, TH1* aEfficiency, TH1*
       if (norm > 0)
       {
         for (Int_t t = kStartBin; t<kMaxT; t++)
-          inverseResponse[t][m] = response[t][m] * prior[t] * efficiency[t] / norm;\r
+          inverseResponse[t][m] = response[t][m] * prior[t] * efficiency[t] / norm;
       }
       else
       {
@@ -725,7 +733,7 @@ Int_t AliUnfolding::UnfoldWithBayesian(TH2* correlation, TH1* aEfficiency, TH1*
   delete[] result;
   delete[] efficiency;
   delete[] binWidths;
-  delete[] normResponse;\r
+  delete[] normResponse;
 
   for (Int_t i=0; i<kMaxT; i++)
   {
@@ -1070,16 +1078,16 @@ Double_t AliUnfolding::RegularizationPowerLaw(TVectorD& params)
     middle = TMath::Power(params[i-1] / fgUnfoldedAxis->GetBinWidth(i),fgPowern);
 
     if(middle>0) {
-      right  = TMath::Power(params[i] / fgUnfoldedAxis->GetBinWidth(i),fgPowern)/middle;
+      right  = TMath::Power(params[i] / fgUnfoldedAxis->GetBinWidth(i+1),fgPowern)/middle;
 
       left   = TMath::Power(params[i-2] / fgUnfoldedAxis->GetBinWidth(i-1),fgPowern)/middle;
       
       middle = 1.;
       
       Double_t der1 = (right - middle) / ((fgUnfoldedAxis->GetBinWidth(i+1) + fgUnfoldedAxis->GetBinWidth(i)) / 2);
-      Double_t der2 = (middle - left) / ((fgUnfoldedAxis->GetBinWidth(i-1) + fgUnfoldedAxis->GetBinWidth(i-2)) / 2);
+      Double_t der2 = (middle - left) / ((fgUnfoldedAxis->GetBinWidth(i) + fgUnfoldedAxis->GetBinWidth(i-1)) / 2);
     
-      chi2 += (der1 - der2) * (der1 - der2)/ ( fgUnfoldedAxis->GetBinWidth(i)/2. + fgUnfoldedAxis->GetBinWidth(i-1) + fgUnfoldedAxis->GetBinWidth(i-2)/2.)/( fgUnfoldedAxis->GetBinWidth(i)/2. + fgUnfoldedAxis->GetBinWidth(i-1) + fgUnfoldedAxis->GetBinWidth(i-2)/2.);// / error;
+      chi2 += (der1 - der2) * (der1 - der2)/ ( fgUnfoldedAxis->GetBinWidth(i+1)/2. + fgUnfoldedAxis->GetBinWidth(i) + fgUnfoldedAxis->GetBinWidth(i-1)/2.)/( fgUnfoldedAxis->GetBinWidth(i+1)/2. + fgUnfoldedAxis->GetBinWidth(i) + fgUnfoldedAxis->GetBinWidth(i-1)/2.);// / error;
       //   printf("i: %d chi2 = %f\n",i,chi2);
     }
 
@@ -1106,14 +1114,13 @@ Double_t AliUnfolding::RegularizationLogLog(TVectorD& params)
     if ((*fgEfficiency)(i-1) == 0 || (*fgEfficiency)(i) == 0 || (*fgEfficiency)(i-2) == 0)
      continue;
 
-
-    Double_t right  = log(params[i] / (*fgEfficiency)(i) / fgUnfoldedAxis->GetBinWidth(i));
-    Double_t middle = log(params[i-1] / (*fgEfficiency)(i-1) / fgUnfoldedAxis->GetBinWidth(i-1));
-    Double_t left   = log(params[i-2] / (*fgEfficiency)(i-2) / fgUnfoldedAxis->GetBinWidth(i-2));
+    Double_t right  = log(params[i] / (*fgEfficiency)(i) / fgUnfoldedAxis->GetBinWidth(i+1));
+    Double_t middle = log(params[i-1] / (*fgEfficiency)(i-1) / fgUnfoldedAxis->GetBinWidth(i));
+    Double_t left   = log(params[i-2] / (*fgEfficiency)(i-2) / fgUnfoldedAxis->GetBinWidth(i-1));
     
     Double_t der1 = (right - middle) / ( log(fgUnfoldedAxis->GetBinCenter(i+1)) - log(fgUnfoldedAxis->GetBinCenter(i)) );
     Double_t der2 = (middle - left) /( log(fgUnfoldedAxis->GetBinCenter(i)) - log(fgUnfoldedAxis->GetBinCenter(i-1)) );
-
+    
     double tmp = (log(fgUnfoldedAxis->GetBinCenter(i+1)) - log(fgUnfoldedAxis->GetBinCenter(i-1)))/2.;
     Double_t dder = (der1-der2) / tmp;
 
@@ -1332,16 +1339,9 @@ void AliUnfolding::DrawResults(TH2* correlation, TH1* efficiency, TH1* measured,
   {
     params[i] = initialConditions->GetBinContent(i+1) * efficiency->GetBinContent(i+1);
 
-    Bool_t fix = kFALSE;
     if (params[i] < 0)
-    {
-      fix = kTRUE;
       params[i] = -params[i];
-    }
     params[i]=TMath::Sqrt(params[i]);
-
-    //cout << "params[" << i << "] " << params[i] << endl;
-
   } 
 
   Double_t chi2;