PWGCF/FEMTOSCOPY/macros/Makec3EAfile.C

   1 #include <math.h>
   2 #include <time.h>
   3 #include <stdio.h>
   4 #include <stdlib.h>
   5 #include <Riostream.h>
   6
   7 #include "TVector2.h"
   8 #include "TFile.h"
   9 #include "TString.h"
  10 #include "TF1.h"
  11 #include "TH1.h"
  12 #include "TH2.h"
  13 #include "TH3.h"
  14 #include "TProfile.h"
  15 #include "TProfile2D.h"
  16 #include "TMath.h"
  17 #include "TText.h"
  18 #include "TRandom3.h"
  19 #include "TArray.h"
  20 #include "TLegend.h"
  21 #include "TStyle.h"
  22 #include "TMinuit.h"
  23 #include "TCanvas.h"
  24 #include "TPad.h"
  25
  26 #define BohrR 1963.6885
  27 #define FmToGeV 0.19733 // conversion to fm
  28 #define PI 3.1415926
  29 #define masspiC 0.1395702 // pi+ mass (GeV/c^2)
  30
  31 using namespace std;
  32
  33
  34
  35 void Makec3EAfile(){
  36
  37
  38   TFile *infile;
  39   TMinuit *fit;
  40   //
  41   TH3D *PbPbEA = new TH3D("PbPbEA","",2,0.5,2.5, 4,0.5,4.5, 50,-0.5,49.5);
  42   PbPbEA->SetDirectory(0);
  43   //
  44   TH3D *pPbEA = new TH3D("pPbEA","",2,0.5,2.5, 4,0.5,4.5, 50,-0.5,49.5);
  45   pPbEA->SetDirectory(0);
  46   //
  47   TH3D *ppEA = new TH3D("ppEA","",2,0.5,2.5, 4,0.5,4.5, 50,-0.5,49.5);
  48   ppEA->SetDirectory(0);
  49   //
  50
  51
  52
  53
  54   //
  55   //////////////////////////////
  56   double value=0, value_e=0;
  57   //
  58   for(int Gindex=0; Gindex<46; Gindex++){
  59     // PbPb
  60     for(int FT=0; FT<2; FT++){// EW or LG
  61
  62       TString *name1 = new TString("FitFiles/FitFile_CT0_FT");
  63       *name1 += FT;
  64       name1->Append("_G");
  65       *name1 += Gindex;
  66       name1->Append(".root");
  67       infile = new TFile(name1->Data(),"READ");
  68       fit = (TMinuit*)infile->Get("MyMinuit_c3");
  69       for(int parNum=0; parNum<4; parNum++){
  70         fit->GetParameter(parNum+1, value,value_e);
  71         PbPbEA->SetBinContent(FT+1, parNum+1, Gindex+1, value);
  72       }
  73       infile->Close();
  74
  75       //
  76       // pPb
  77       TString *name1 = new TString("FitFiles/FitFile_CT1_FT");
  78       *name1 += FT;
  79       name1->Append("_G");
  80       *name1 += Gindex;
  81       name1->Append(".root");
  82       infile = new TFile(name1->Data(),"READ");
  83       fit = (TMinuit*)infile->Get("MyMinuit_c3");
  84       for(int parNum=0; parNum<4; parNum++){
  85         fit->GetParameter(parNum+1, value,value_e);
  86         pPbEA->SetBinContent(FT+1, parNum+1, Gindex+1, value);
  87       }
  88       infile->Close();
  89       //
  90       //
  91       TString *name1 = new TString("FitFiles/FitFile_CT2_FT");
  92       *name1 += FT;
  93       name1->Append("_G");
  94       *name1 += Gindex;
  95       name1->Append(".root");
  96       infile = new TFile(name1->Data(),"READ");
  97       fit = (TMinuit*)infile->Get("MyMinuit_c3");
  98       for(int parNum=0; parNum<4; parNum++){
  99         fit->GetParameter(parNum+1, value,value_e);
 100         ppEA->SetBinContent(FT+1, parNum+1, Gindex+1, value);
 101       }
 102       infile->Close();
 103     }
 104   }
 105   // blank for the rest
 106   for(int Gindex=46; Gindex<50; Gindex++){
 107     for(int FT=0; FT<2; FT++){// EW or LG
 108       for(int parNum=0; parNum<4; parNum++){
 109         PbPbEA->SetBinContent(FT+1, parNum+1, Gindex+1, PbPbEA->GetBinContent(FT+1, parNum+1, 46));
 110         pPbEA->SetBinContent(FT+1, parNum+1, Gindex+1, pPbEA->GetBinContent(FT+1, parNum+1, 46));
 111         ppEA->SetBinContent(FT+1, parNum+1, Gindex+1, ppEA->GetBinContent(FT+1, parNum+1, 46));
 112       }
 113     }
 114   }
 115
 116   // Convert Lam_3 to proper EA normalization
 117   for(int Gindex=0; Gindex<50; Gindex++){
 118     for(int FT=0; FT<2; FT++){// EW or LG
 119       PbPbEA->SetBinContent(FT+1, 1, Gindex+1, pow(PbPbEA->GetBinContent(FT+1, 1, Gindex+1)/ 2., 1/3.));
 120       pPbEA->SetBinContent(FT+1, 1, Gindex+1, pow(pPbEA->GetBinContent(FT+1, 1, Gindex+1)/ 2., 1/3.));
 121       ppEA->SetBinContent(FT+1, 1, Gindex+1, pow(ppEA->GetBinContent(FT+1, 1, Gindex+1)/ 2., 1/3.));
 122     }
 123   }
 124
 125
 126   TFile *outfile=new TFile("c3EAfile_temp.root","RECREATE");
 127   PbPbEA->Write();
 128   pPbEA->Write();
 129   ppEA->Write();
 130   //
 131   outfile->Close();
 132
 133
 134 }