Adding the new histograms in merging of Event cuts and new way of nomralization in...
authormchojnac <mchojnac@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 6 May 2013 09:16:00 +0000 (09:16 +0000)
committermchojnac <mchojnac@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 6 May 2013 09:16:00 +0000 (09:16 +0000)
PWGLF/SPECTRA/PiKaPr/TestAOD/AliSpectraBothEventCuts.cxx
PWGLF/SPECTRA/PiKaPr/TestAOD/AnalysisBoth.C
PWGLF/SPECTRA/PiKaPr/TestAOD/QAPlotsBoth.C

index 10703f8..f3f316c 100644 (file)
@@ -436,6 +436,9 @@ Long64_t AliSpectraBothEventCuts::Merge(TCollection* list)
   // TList collections_histoQVectorNeg;
   TList collections_histoQVector;
   TList collections_histoEP;
+  TList collections_histoVtxAftSelwithoutZvertexCut;
+  TList collections_histoVtxalltriggerEventswithMCz;
+  TList collections_histoVtxAftSelwithoutZvertexCutusingMCz;                   
 
   Int_t count = 0;
 
@@ -464,6 +467,13 @@ Long64_t AliSpectraBothEventCuts::Merge(TCollection* list)
     collections_histoQVector.Add(histo_histoQVector);
     TH1F * histo_histoEP = entry->GetHistoEP();      
     collections_histoEP.Add(histo_histoEP);
+    TH1F* histo_histoVtxAftSelwithoutZvertexCut=entry->GetHistoVtxAftSelwithoutZvertexCut();
+    collections_histoVtxAftSelwithoutZvertexCut.Add(histo_histoVtxAftSelwithoutZvertexCut);
+    TH1F* histo_histoVtxalltriggerEventswithMCz=entry->GetHistoVtxGenerated();
+    collections_histoVtxalltriggerEventswithMCz.Add(histo_histoVtxalltriggerEventswithMCz);
+    TH1F* histo_histoVtxAftSelwithoutZvertexCutusingMCz=entry->GetHistoVtxAftSelwithoutZvertexCutusingMCz();
+    collections_histoVtxAftSelwithoutZvertexCutusingMCz.Add(histo_histoVtxAftSelwithoutZvertexCutusingMCz);    
+
     count++;
   }
   
@@ -477,7 +487,11 @@ Long64_t AliSpectraBothEventCuts::Merge(TCollection* list)
   // fHistoQVectorNeg->Merge(&collections_histoQVectorNeg);
   fHistoQVector->Merge(&collections_histoQVector);
   fHistoEP->Merge(&collections_histoEP);
-  
+
+  fHistoVtxAftSelwithoutZvertexCut->Merge(&collections_histoVtxAftSelwithoutZvertexCut);
+  fHistoVtxalltriggerEventswithMCz->Merge(&collections_histoVtxalltriggerEventswithMCz);
+  fHistoVtxAftSelwithoutZvertexCutusingMCz->Merge(&collections_histoVtxAftSelwithoutZvertexCutusingMCz);
 
   delete iter;
 
index 988e78c..3906b47 100644 (file)
@@ -30,7 +30,10 @@ enum {
  kveretxcorrectionandbadchunkscorr=0x10,
  kmcisusedasdata=0x20,
  kdonotusedcacuts=0x40,
- kuseprimaryPIDcont=0x80                       
+ kuseprimaryPIDcont=0x80,
+ knormalizationwithbin0integralsdata=0x100,
+ knormalizationwithbin0integralsMC=0x200,
+                                       
 };     
 
 Bool_t OpenFile(TString dirname, TString outputname, Bool_t mcflag,Bool_t mcasdata=false);
@@ -119,24 +122,34 @@ void AnalysisBoth (UInt_t options=0xF,TString outdate, TString outnamedata, TStr
        GetPtHistFromPtDCAhisto("hHistPtRecSigmaPrimary","conPIDprimary",managermc,contPIDpri,dcacutxy);
 
        
-//     Double_t neventsdata =  ecutsdata->NumberOfEvents();
        Double_t neventsmcall = 1 ;  //if loop over MC is done after or befor events cuts this will be changed 
        Double_t neventsdata =  1;
        Double_t neventsmc =  1;
 
+       //Normaliztion of MCtruth depends if the loop was done after of before ESD event cuts.
+       //In currect code this cannot be check on the level of macro.
+       //If the loop was done before MC should be done to all processed events (NumberOfProcessedEvents())
+       //If loop was done after MC should be normalized to all accepted events (NumberOfEvents()) 
+       // The option one will be alaways use.
+       
+       neventsmcall= ecutsmc->NumberOfProcessedEvents();
        if(options&knormalizationtoeventspassingPhySel)
        {
-               neventsmcall= ecutsmc->NumberOfProcessedEvents();
+               //neventsmcall= ecutsmc->NumberOfProcessedEvents();
                 neventsdata=ecutsdata->NumberOfPhysSelEvents();
                 neventsmc=ecutsmc->NumberOfPhysSelEvents();
        }
+       else if ((options&knormalizationwithbin0integralsdata)||(options&knormalizationwithbin0integralsMC))
+       {
+               neventsdata=Normaliztionwithbin0integrals(options);
+               neventsmc=ecutsmc->NumberOfPhysSelEvents();
+       }
        else
        {
-               neventsdata=ecutsdata->NumberOfEvents();
+               //neventsdata=ecutsdata->NumberOfEvents(); //number of accepted events
                 neventsmc=ecutsmc->NumberOfEvents();
                neventsmcall= ecutsmc->NumberOfEvents();
 
-
        }
        GetMCTruth(MCTruth);
        
@@ -314,7 +327,11 @@ void AnalysisBoth (UInt_t options=0xF,TString outdate, TString outnamedata, TStr
        TList* listqa=new TList();
        TList* canvaslist=new TList();
        Float_t vertexcorrection=1.0;
-       Float_t corrbadchunksvtx=QAPlotsBoth(managerdata,managermc,ecutsdata,ecutsmc,tcutsdata,tcutsmc,listqa,canvaslist);
+       Float_t corrbadchunksvtx=1.0;
+       if ((options&knormalizationwithbin0integralsdata)||(options&knormalizationwithbin0integralsMC))
+               corrbadchunksvtx=QAPlotsBoth(managerdata,managermc,ecutsdata,ecutsmc,tcutsdata,tcutsmc,listqa,canvaslist,0);
+       else
+               corrbadchunksvtx=QAPlotsBoth(managerdata,managermc,ecutsdata,ecutsmc,tcutsdata,tcutsmc,listqa,canvaslist,1);
        if (options&kveretxcorrectionandbadchunkscorr)
                vertexcorrection=corrbadchunksvtx;
        cout<<" VTX corr="<<vertexcorrection<<endl;
@@ -341,6 +358,7 @@ void AnalysisBoth (UInt_t options=0xF,TString outdate, TString outnamedata, TStr
        canvaslist->Write("outputcanvas",TObject::kSingleKey);
 
        fout->Close();
+       //Normaliztionwithbin0integrals();
 
 }
 
@@ -1168,4 +1186,60 @@ Short_t DCAfitsettings (Float_t pt, Int_t type)
                value=value+1;
        return value;   
 
-}  
+} 
+
+Float_t Normaliztionwithbin0integrals(UInt_t options)
+{
+       
+       TH1F* bin0mcRec=(TH1F*)ecutsmc->GetHistoVtxGenerated()->Clone("Bin0_rec");
+       TH1F* bin0mcMC=(TH1F*)ecutsmc->GetHistoVtxGenerated()->Clone("Bin0_MC");
+
+       TH1F* vertexmc=ecutsmc->GetHistoVtxAftSelwithoutZvertexCut(); 
+       TH1F* vertexmcMCz=ecutsmc->GetHistoVtxAftSelwithoutZvertexCutusingMCz(); 
+       TH1F* vertexdata=ecutsdata->GetHistoVtxAftSelwithoutZvertexCut();
+
+       TH1I* histodata=ecutsdata->GetHistoCuts();
+       TH1I* histomc=ecutsmc->GetHistoCuts();
+
+       Float_t dataevents=(Float_t)histodata->GetBinContent(3);
+       cout<<histodata->GetBinContent(2)<<endl;
+       Float_t databin0events=((Float_t)histodata->GetBinContent(2))-((Float_t)histodata->GetBinContent(4));   
+
+       bin0mcRec->Sumw2();
+       bin0mcMC->Sumw2();
+               
+       bin0mcRec->Add(vertexmc,-1);
+       bin0mcMC->Add(vertexmcMCz,-1);
+       
+       bin0mcRec->Divide(vertexmc);
+       bin0mcMC->Divide(vertexmcMCz);
+       
+       bin0mcRec->Multiply(vertexdata);
+       bin0mcMC->Multiply(vertexdata);
+       
+       Float_t bin0mcRecN=0.0;
+       Float_t bin0mcMCN=0.0;
+
+       for (int i=0;i<=bin0mcRec->GetXaxis()->GetNbins();i++)
+       {
+               bin0mcRecN+=bin0mcRec->GetBinContent(i);
+               bin0mcMCN+=bin0mcMC->GetBinContent(i);
+
+       }
+       bin0mcRec->Scale(databin0events/bin0mcRecN);
+       bin0mcMC->Scale(databin0events/bin0mcMCN);              
+       
+       Int_t binmin=bin0mcRec->GetXaxis()->FindBin(-10);
+       Int_t binmax=bin0mcRec->GetXaxis()->FindBin(10)-1;
+       cout<<  bin0mcRec->GetXaxis()->GetBinLowEdge(binmin)<<" "<<bin0mcRec->GetXaxis()->GetBinUpEdge(binmax)<<endl;
+       cout<<bin0mcRecN<<" "<<bin0mcMCN<<" "<<databin0events<<endl;    
+       cout<<dataevents<<" normalization "<<dataevents+bin0mcRec->Integral(binmin,binmax)<<" "<<dataevents+bin0mcMC->Integral(binmin,binmax)<<endl;
+       cout<<histodata->GetBinContent(2)<<" "<<histodata->GetBinContent(4)<<endl;
+       if ((options&knormalizationwithbin0integralsdata)==knormalizationwithbin0integralsdata)
+               return  dataevents+bin0mcRec->Integral(binmin,binmax);
+       else if ((options&kknormalizationwithbin0integralsMC)==knormalizationwithbin0integralsdata)
+               return dataevents+bin0mcMC->Integral(binmin,binmax) ;
+       else
+               return 1;               
+}
index df1f313..9ff796a 100644 (file)
@@ -1,7 +1,7 @@
 Float_t QAPlotsBoth( AliSpectraBothHistoManager* hman_data, AliSpectraBothHistoManager* hman_mc,
              AliSpectraBothEventCuts* ecuts_data, AliSpectraBothEventCuts* ecuts_mc,
              AliSpectraBothTrackCuts* tcuts_data, AliSpectraBothTrackCuts* tcuts_mc,
-             TList * flistqa,TList * flistcanvas)
+             TList * flistqa,TList * flistcanvas,Bool_t fullicorr=kTRUE)
 {
 TString pidmethods[3]={"TPC","TOF","TPCTOF"};  
        Double_t neventsdata =  ecutsdata->NumberOfPhysSelEvents();
@@ -96,6 +96,8 @@ TString pidmethods[3]={"TPC","TOF","TPCTOF"};
                 binstartx++;
         binstartx++;
        flistcanvas->Add(plot_on_canvas("vertex",fHistoVtxAftSeldata,fHistoVtxAftSelmc));
+       /*
+
        TF1* fdata=new TF1("dataveretxfit","gausn");
        TF1* fmc=new TF1("mcveretxfit","gausn");
        //we strat fit a second not empty bin
@@ -105,8 +107,20 @@ TString pidmethods[3]={"TPC","TOF","TPCTOF"};
        fHistoVtxAftSelmc->Fit("mcveretxfit","0","",minfit,-1.0*minfit);
        Float_t datavertexratio=fHistoVtxAftSeldata->Integral(-1,-1,"width")/fdata->GetParameter(0);
        Float_t mcvertexratio=fHistoVtxAftSelmc->Integral(-1,-1,"width")/fmc->GetParameter(0);
+       */
+       //Event cut histo 
+       TH1I* histodata=ecuts_data->GetHistoCuts();
+       TH1I* histomc=ecuts_mc->GetHistoCuts();
        
+       Int_t events_data=histodata->GetBinContent(3);
+       Int_t events_mc=histomc->GetBinContent(3);
+
+       if(events_data==0&&events_mc==0)
+               return 0;
 
+
+       Float_t datavertexratio=((Float_t)(events_data))/((Float_t)histodata->GetBinContent(4));
+         Float_t mcvertexratio=((Float_t)(events_mc))/((Float_t)histomc->GetBinContent(4));
         TH1F* fHistoEtaAftSeldata=(TH1F*)ecuts_data->GetHistoEtaAftSel();
         TH1F* fHistoEtaAftSelmc=(TH1F*)ecuts_mc->GetHistoEtaAftSel();
        flistcanvas->Add(plot_on_canvas("ETA",fHistoEtaAftSeldata,fHistoEtaAftSelmc));
@@ -140,9 +154,15 @@ TString pidmethods[3]={"TPC","TOF","TPCTOF"};
        cout<<"Bad chunks "<<badchunksfraction<<endl;   
        binzero->Draw("E1");
        flistcanvas->Add(cbc);
+       if(TMath::Abs(hmul->GetEntries()/events_mc-1.0)>0.001)
+               cout<<"MC merging problem"<<endl;
+       if(TMath::Abs(hman_data->GetGenMulvsRawMulHistogram("hHistGenMulvsRawMul")->GetEntries()/events_data-1.0)>0.001)
+               cout<<"Data merging problem"<<endl;
        
-       return (1.0-badchunksfraction)*mcvertexratio/datavertexratio;
-
+       if(fullicorr)
+               return (1.0-badchunksfraction)*mcvertexratio/datavertexratio;
+       else 
+               return (1.0-badchunksfraction)*mcvertexratio;
 }
 TCanvas* plot_on_canvas(TString name, TH1* h1,TH1* h2)
 {