2D case

CMakeLists.txt

@@ -8,6 +8,7 @@ set(SOURCES
         src/BordersFinder.cpp
         src/BordersFinderHelper.cpp
         src/Getter.cpp
+        src/Converter1D.cpp
         glauber/Fitter.cpp 
         glauber/FitterHelper.cpp 
         )

macro/TestGetter.C

@@ -2,7 +2,7 @@
 void TestGetter(const TString file = "../build/test.root", const TString infile = "../input/test_input.root" )
 {
     std::unique_ptr <TFile> f {TFile::Open(file, "read")};    
-    std::unique_ptr <Centrality::Getter> getter {(Centrality::Getter*)f->Get("getter")};
+    std::unique_ptr <Centrality::Getter> getter {(Centrality::Getter*)f->Get("centr_getter")};
     
     const float value = 20.;
     const float centrality = getter->GetCentrality(value);

main.cpp

 #include "Getter.h"
 #include "BordersFinder.h"
+#include "Converter1D.h"
 
 #include <iostream>
 #include <vector>
@@ -12,6 +13,7 @@
 #include <TROOT.h>
 #include "TTree.h"
 #include "TFile.h"
+#include "TH2.h"
 
 
 int main(int argc, char **argv) {
@@ -20,9 +22,19 @@ int main(int argc, char **argv) {
 
   std::unique_ptr <TFile> fIn {TFile::Open(argv[1], "read")};
   std::unique_ptr <TH1F> histo {(TH1F*) (fIn->Get(argv[2]))};
-    
+
+  std::unique_ptr <TH2F> histo2d {(TH2F*) (fIn->Get("hMEcorr"))};
+  
+  Centrality::Converter1D conv2d;
+  conv2d.SetHisto2D( std::move(*histo2d) );
+  conv2d.Convert();
+  
+
+  if (false)
+  {
+  
   Centrality::BordersFinder bf;
-  bf.SetHisto(*histo);
+  bf.SetHisto( *histo);
   bf.SetRanges( 10,0,100 );   // number of bins, min, max value
 //   bf.SetRanges( {0,10,30,60,100} );  // centrality bins borders with array
   bf.IsSpectator(true);  // true if impact parameter b correlated with estimator (spectators eneggy), false - anticorrelated (multiplicity of produced particles) 
@@ -31,7 +43,7 @@ int main(int argc, char **argv) {
   
   std::string outfilename = "test.root";
   bf.SaveBorders(outfilename);
-  
+  }
   auto end = std::chrono::system_clock::now();
   std::chrono::duration<double> elapsed_seconds = end - start;
   std::cout << "elapsed time: " << elapsed_seconds.count() << " s\n";

src/BordersFinder.cpp

@@ -52,7 +52,7 @@ void BordersFinder::SaveBorders(std::string filename)
     h.QA(getter, histo_);
     h.PlotHisto(getter, histo_);
     
-    getter.Write("getter");
+    getter.Write("centr_getter");
     
     f->Close();
 }

src/Converter1D.cpp

+#include "Converter1D.h"
+#include <iostream>
+
+#include "TProfile.h"
+#include "TMath.h"
+
+
+// ClassImp(Centrality::Converter1D)
+
+namespace Centrality {
+
+void Converter1D::Convert()
+{
+    std::cout << histo2d_.GetNbinsX() << std::endl;
+    std::cout << histo2d_.GetNbinsY() << std::endl;
+}
+    
+    
+std::unique_ptr<TF1> Converter1D::Fit2D()
+{
+    std::unique_ptr<TProfile> prof{ histo2d_.ProfileX() };
+    std::unique_ptr<TF1> fitf = new TF1("fit", "pol2", histo2d_.GetXaxis()->GetXmin(), histo2d_.GetXaxis()->GetXmax() );
+    
+    prof->Fit(fitf);
+    return fitf;
+}
+    
+    
+/**
+ *
+ * @param par parameters of polinomial function
+ * @param x argument
+ * @param kb return value
+ * @param N order
+ */
+void Converter1D::FindNorm (float par[], float x, float kb[], int N)
+{
+    const float dx = (histo2d_.GetXaxis()->GetXmax() - histo2d_.GetXaxis()->GetXmin()) / 10000. ;
+
+    /* left */
+    const float y1 = polN(par, x - dx, N);
+    /* right */
+    const float y2 = polN(par, x + dx, N);
+    
+    // cx*x + cy*y + c == 0
+
+    /* 1/df */
+    const float cx = 1/(y2 - y1);
+    /* 1/2dx */
+    const float cy = 0.5/dx;
+        
+    const float c = -cx*x - cy*polN(par, x, N);
+    
+    kb[0] = -cx/cy;
+    kb[1] = -c/cy;
+}    
+    
+
+}

src/Converter1D.h

+//
+// Created by Viktor Klochkov on 08.08.18.
+//
+
+#ifndef CENTRALITY_CONVERTER1D_H
+#define CENTRALITY_CONVERTER1D_H
+
+#include "vector"
+
+#include "TH2.h"
+#include "TH1.h"
+
+
+namespace Centrality {
+    
+class Converter1D {
+
+public:
+    
+    Converter1D(){}
+        
+    void SetHisto2D(TH2F&& histo2d) { histo2d_ = histo2d; }
+    TH2F&& GetHisto2D() { return std::move(histo2d_); }
+    
+    void Convert();
+    std::unique_ptr<TF1> Fit2D();
+
+    void FindNorm (float par[], float x, float kb[], int N);
+
+
+    /**
+     *
+     * @param par parameters of polinom
+     * @param x argument
+     * @param N order
+     * @return
+     */
+    float polN (float par[], float x, uint N)
+    {
+        float res = par[0];
+        float xn = 1.;
+        for (uint i=1; i<N;  ++i){
+            xn *= x;
+            res += par[i]*xn;
+        }
+        return res;
+    }
+
+private:
+    
+    TH2F histo2d_;
+    TH1F histo1d_;
+    
+    
+//     ClassDef(Converter1D, 1);
+    
+};
+}
+
+#endif //CENTRALITY_CONVERTER1D_H

src/Getter.h

@@ -41,7 +41,6 @@ private:
     bool isspectator_{false};
 
     ClassDef(Getter, 2);
-
     
 };
 }