ANALYSIS/AliPPVsMultUtils.cxx

   1 /*****************************************************************
   2
   3  General utility class for doing pp vs multiplicity studies
   4
   5  --- Functionality being added on demand.
   6
   7  --- Please report any bugs, comments, suggestions to:
   8  david.dobrigkeit.chinellato@cern.ch
   9
  10  *****************************************************************/
  11
  12 #include "AliVEvent.h"
  13 #include "AliESDEvent.h"
  14 #include "AliAODEvent.h"
  15 #include "AliVVertex.h"
  16 #include "AliLog.h"
  17 #include "AliAODVertex.h"
  18 #include "AliVTrack.h"
  19 #include "AliVEvent.h"
  20 #include <TMatrixDSym.h>
  21 #include <TMath.h>
  22 #include "AliESDUtils.h"
  23 #include "AliPPVsMultUtils.h"
  24 #include <TFile.h>
  25
  26
  27 ClassImp(AliPPVsMultUtils)
  28
  29 //______________________________________________________________________
  30 AliPPVsMultUtils::AliPPVsMultUtils():TObject(),
  31 fRunNumber(0),
  32 fBoundaryHisto_V0M(0),
  33 fBoundaryHisto_V0A(0),
  34 fBoundaryHisto_V0C(0),
  35 fBoundaryHisto_V0MEq(0),
  36 fBoundaryHisto_V0AEq(0),
  37 fBoundaryHisto_V0CEq(0)
  38 {
  39     // Default contructor
  40 }
  41
  42
  43 //_____________________________________________________________________________
  44 AliPPVsMultUtils::AliPPVsMultUtils(const AliPPVsMultUtils &c) : TObject(c),
  45 fRunNumber(0),
  46 fBoundaryHisto_V0M(0),
  47 fBoundaryHisto_V0A(0),
  48 fBoundaryHisto_V0C(0),
  49 fBoundaryHisto_V0MEq(0),
  50 fBoundaryHisto_V0AEq(0),
  51 fBoundaryHisto_V0CEq(0)
  52 {
  53     //
  54     // copy constructor - untested
  55     //
  56     ((AliPPVsMultUtils &) c).Copy(*this);
  57 }
  58
  59 //_____________________________________________________________________________
  60 AliPPVsMultUtils &AliPPVsMultUtils::operator=(const AliPPVsMultUtils &c)
  61 {
  62     //
  63     // Assignment operator - untested
  64     //
  65
  66     if (this != &c) ((AliPPVsMultUtils &) c).Copy(*this);
  67     return *this;
  68 }
  69
  70 //______________________________________________________________________
  71 Float_t AliPPVsMultUtils::GetMultiplicityPercentile(AliESDEvent *event, TString lMethod)
  72 // Function to get multiplicity quantiles in ESDs. All methods would work
  73 // with AliVEvent except GetPrimaryVertexSPD which exists in both ESD and AOD
  74 // but not in the AliVEvent class; therefore two functions are included.
  75 // N.B.: Any change to this method has to be propagated by hand to the function
  76 // below!
  77 {
  78     Int_t lRequestedRunNumber = event->GetRunNumber();
  79     if( lRequestedRunNumber != fRunNumber ) LoadCalibration( lRequestedRunNumber );
  80
  81     Float_t lreturnval = -1;
  82
  83     //Get VZERO Information for multiplicity later
  84     AliVVZERO* esdV0 = event->GetVZEROData();
  85     if (!esdV0) {
  86         AliError("AliVVZERO not available");
  87         return -1;
  88     }
  89
  90     // VZERO PART
  91     Float_t  multV0A  = 0;            //  multiplicity from V0 reco side A
  92     Float_t  multV0C  = 0;            //  multiplicity from V0 reco side C
  93     Float_t  multV0AEq  = 0;          //  multiplicity from V0 reco side A
  94     Float_t  multV0CEq  = 0;          //  multiplicity from V0 reco side C
  95     Float_t  multV0ACorr  = 0;            //  multiplicity from V0 reco side A
  96     Float_t  multV0CCorr  = 0;            //  multiplicity from V0 reco side C
  97
  98     //Non-Equalized Signal: copy of multV0ACorr and multV0CCorr from AliCentralitySelectionTask
  99     //Getters for uncorrected multiplicity
 100     multV0A=esdV0->GetMTotV0A();
 101     multV0C=esdV0->GetMTotV0C();
 102     const AliESDVertex *lPrimarySPDVtx         = event->GetPrimaryVertexSPD();
 103
 104     //Get Z vertex position of SPD vertex (why not Tracking if available?)
 105     Float_t zvtx = lPrimarySPDVtx->GetZ();
 106
 107     //Acquire Corrected multV0A
 108     multV0ACorr = AliESDUtils::GetCorrV0A(multV0A,zvtx);
 109     multV0CCorr = AliESDUtils::GetCorrV0C(multV0C,zvtx);
 110
 111     // Equalized signals // From AliCentralitySelectionTask // Updated
 112     for(Int_t iCh = 32; iCh < 64; ++iCh) {
 113         Double_t mult = event->GetVZEROEqMultiplicity(iCh);
 114         multV0AEq += mult;
 115     }
 116     for(Int_t iCh = 0; iCh < 32; ++iCh) {
 117         Double_t mult = event->GetVZEROEqMultiplicity(iCh);
 118         multV0CEq += mult;
 119     }
 120
 121     if ( lMethod == "V0M" ) lreturnval =
 122         fBoundaryHisto_V0M -> GetBinContent( fBoundaryHisto_V0M->FindBin(multV0ACorr+multV0CCorr) );
 123     if ( lMethod == "V0A" ) lreturnval =
 124         fBoundaryHisto_V0A -> GetBinContent( fBoundaryHisto_V0A->FindBin(multV0ACorr) );
 125     if ( lMethod == "V0C" ) lreturnval =
 126         fBoundaryHisto_V0C -> GetBinContent( fBoundaryHisto_V0C->FindBin(multV0CCorr) );
 127     //equalized
 128     if ( lMethod == "V0MEq" ) lreturnval =
 129         fBoundaryHisto_V0MEq -> GetBinContent( fBoundaryHisto_V0MEq->FindBin(multV0AEq+multV0CEq) );
 130     if ( lMethod == "V0AEq" ) lreturnval =
 131         fBoundaryHisto_V0AEq -> GetBinContent( fBoundaryHisto_V0AEq->FindBin(multV0AEq) );
 132     if ( lMethod == "V0CEq" ) lreturnval =
 133         fBoundaryHisto_V0CEq -> GetBinContent( fBoundaryHisto_V0CEq->FindBin(multV0CEq) );
 134
 135     return lreturnval;
 136 }
 137
 138 //______________________________________________________________________
 139 Float_t AliPPVsMultUtils::GetMultiplicityPercentile(AliAODEvent *event, TString lMethod)
 140 // Carbon-copy version of ESD function (GetPrimaryVertexSPD won't work
 141 // in the base class AliVEvent...FIXME)
 142 {
 143     Int_t lRequestedRunNumber = event->GetRunNumber();
 144     if( lRequestedRunNumber != fRunNumber ) LoadCalibration( lRequestedRunNumber );
 145
 146     Float_t lreturnval = -1;
 147
 148     //Get VZERO Information for multiplicity later
 149     AliVVZERO* esdV0 = event->GetVZEROData();
 150     if (!esdV0) {
 151         AliError("AliVVZERO not available");
 152         return -1;
 153     }
 154
 155     // VZERO PART
 156     Float_t  multV0A  = 0;            //  multiplicity from V0 reco side A
 157     Float_t  multV0C  = 0;            //  multiplicity from V0 reco side C
 158     Float_t  multV0AEq  = 0;          //  multiplicity from V0 reco side A
 159     Float_t  multV0CEq  = 0;          //  multiplicity from V0 reco side C
 160     Float_t  multV0ACorr  = 0;            //  multiplicity from V0 reco side A
 161     Float_t  multV0CCorr  = 0;            //  multiplicity from V0 reco side C
 162
 163     //Non-Equalized Signal: copy of multV0ACorr and multV0CCorr from AliCentralitySelectionTask
 164     //Getters for uncorrected multiplicity
 165     multV0A=esdV0->GetMTotV0A();
 166     multV0C=esdV0->GetMTotV0C();
 167     const AliAODVertex *lPrimarySPDVtx         = event->GetPrimaryVertexSPD();
 168
 169     //Get Z vertex position of SPD vertex (why not Tracking if available?)
 170     Float_t zvtx = lPrimarySPDVtx->GetZ();
 171
 172     //Acquire Corrected multV0A
 173     multV0ACorr = AliESDUtils::GetCorrV0A(multV0A,zvtx);
 174     multV0CCorr = AliESDUtils::GetCorrV0C(multV0C,zvtx);
 175
 176     // Equalized signals // From AliCentralitySelectionTask // Updated
 177     for(Int_t iCh = 32; iCh < 64; ++iCh) {
 178         Double_t mult = event->GetVZEROEqMultiplicity(iCh);
 179         multV0AEq += mult;
 180     }
 181     for(Int_t iCh = 0; iCh < 32; ++iCh) {
 182         Double_t mult = event->GetVZEROEqMultiplicity(iCh);
 183         multV0CEq += mult;
 184     }
 185
 186     if ( lMethod == "V0M" ) lreturnval =
 187         fBoundaryHisto_V0M -> GetBinContent( fBoundaryHisto_V0M->FindBin(multV0ACorr+multV0CCorr) );
 188     if ( lMethod == "V0A" ) lreturnval =
 189         fBoundaryHisto_V0A -> GetBinContent( fBoundaryHisto_V0A->FindBin(multV0ACorr) );
 190     if ( lMethod == "V0C" ) lreturnval =
 191         fBoundaryHisto_V0C -> GetBinContent( fBoundaryHisto_V0C->FindBin(multV0CCorr) );
 192     //equalized
 193     if ( lMethod == "V0MEq" ) lreturnval =
 194         fBoundaryHisto_V0MEq -> GetBinContent( fBoundaryHisto_V0MEq->FindBin(multV0AEq+multV0CEq) );
 195     if ( lMethod == "V0AEq" ) lreturnval =
 196         fBoundaryHisto_V0AEq -> GetBinContent( fBoundaryHisto_V0AEq->FindBin(multV0AEq) );
 197     if ( lMethod == "V0CEq" ) lreturnval =
 198         fBoundaryHisto_V0CEq -> GetBinContent( fBoundaryHisto_V0CEq->FindBin(multV0CEq) );
 199
 200     return lreturnval;
 201 }
 202
 203 //______________________________________________________________________
 204 Bool_t AliPPVsMultUtils::LoadCalibration(Int_t lLoadThisCalibration)
 205 //To be called if starting analysis on a new run
 206 {
 207     //If Histograms exist, de-allocate to prevent memory leakage
 208     if( fBoundaryHisto_V0M ) {fBoundaryHisto_V0M->Delete(); fBoundaryHisto_V0M = 0x0; }
 209     if( fBoundaryHisto_V0A ) {fBoundaryHisto_V0A->Delete(); fBoundaryHisto_V0A = 0x0; }
 210     if( fBoundaryHisto_V0C ) {fBoundaryHisto_V0C->Delete(); fBoundaryHisto_V0C = 0x0; }
 211     if( fBoundaryHisto_V0MEq ) {fBoundaryHisto_V0MEq->Delete(); fBoundaryHisto_V0MEq = 0x0; }
 212     if( fBoundaryHisto_V0AEq ) {fBoundaryHisto_V0AEq->Delete(); fBoundaryHisto_V0AEq = 0x0; }
 213     if( fBoundaryHisto_V0CEq ) {fBoundaryHisto_V0CEq->Delete(); fBoundaryHisto_V0CEq = 0x0; }
 214
 215     AliInfo(Form( "Loading calibration file for run %i",lLoadThisCalibration) );
 216     TFile *lCalibFile_V0M = 0x0;
 217     TFile *lCalibFile_V0A = 0x0;
 218     TFile *lCalibFile_V0C = 0x0;
 219     TFile *lCalibFile_V0MEq = 0x0;
 220     TFile *lCalibFile_V0AEq = 0x0;
 221     TFile *lCalibFile_V0CEq = 0x0;
 222
 223     //AliInfo("Calling TFile::Open");
 224     lCalibFile_V0M = TFile::Open("$ALICE_ROOT/PWGLF/STRANGENESS/Cascades/corrections/calibration_V0M.root");
 225     lCalibFile_V0A = TFile::Open("$ALICE_ROOT/PWGLF/STRANGENESS/Cascades/corrections/calibration_V0A.root");
 226     lCalibFile_V0C = TFile::Open("$ALICE_ROOT/PWGLF/STRANGENESS/Cascades/corrections/calibration_V0C.root");
 227     lCalibFile_V0MEq = TFile::Open("$ALICE_ROOT/PWGLF/STRANGENESS/Cascades/corrections/calibration_V0MEq.root");
 228     lCalibFile_V0AEq = TFile::Open("$ALICE_ROOT/PWGLF/STRANGENESS/Cascades/corrections/calibration_V0AEq.root");
 229     lCalibFile_V0CEq = TFile::Open("$ALICE_ROOT/PWGLF/STRANGENESS/Cascades/corrections/calibration_V0CEq.root");
 230
 231     //AliInfo("Casting");
 232     //check memory consumption later...
 233     fBoundaryHisto_V0M   = dynamic_cast<TH1F *>(lCalibFile_V0M  -> Get(Form("histocalib%i",lLoadThisCalibration)) );
 234     fBoundaryHisto_V0A   = dynamic_cast<TH1F *>(lCalibFile_V0A  -> Get(Form("histocalib%i",lLoadThisCalibration)) );
 235     fBoundaryHisto_V0C   = dynamic_cast<TH1F *>(lCalibFile_V0C  -> Get(Form("histocalib%i",lLoadThisCalibration)) );
 236     fBoundaryHisto_V0MEq   = dynamic_cast<TH1F *>(lCalibFile_V0MEq  -> Get(Form("histocalib%i",lLoadThisCalibration)) );
 237     fBoundaryHisto_V0AEq   = dynamic_cast<TH1F *>(lCalibFile_V0AEq  -> Get(Form("histocalib%i",lLoadThisCalibration)) );
 238     fBoundaryHisto_V0CEq   = dynamic_cast<TH1F *>(lCalibFile_V0CEq  -> Get(Form("histocalib%i",lLoadThisCalibration)) );
 239
 240     if ( !fBoundaryHisto_V0M   || !fBoundaryHisto_V0A   || !fBoundaryHisto_V0C ||
 241         !fBoundaryHisto_V0MEq || !fBoundaryHisto_V0AEq || !fBoundaryHisto_V0CEq ){
 242         AliFatal(Form("No calibration for run %i exists at the moment!",lLoadThisCalibration));
 243         return kFALSE; //return denial
 244     }
 245
 246     //Careful with manual cleanup if needed: to be implemented
 247     fBoundaryHisto_V0M->SetDirectory(0);
 248     fBoundaryHisto_V0A->SetDirectory(0);
 249     fBoundaryHisto_V0C->SetDirectory(0);
 250     fBoundaryHisto_V0MEq->SetDirectory(0);
 251     fBoundaryHisto_V0AEq->SetDirectory(0);
 252     fBoundaryHisto_V0CEq->SetDirectory(0);
 253
 254     //AliInfo("Closing");
 255
 256     if( lCalibFile_V0M ) lCalibFile_V0M->Close();
 257     if( lCalibFile_V0A ) lCalibFile_V0A->Close();
 258     if( lCalibFile_V0C ) lCalibFile_V0C->Close();
 259     if( lCalibFile_V0MEq ) lCalibFile_V0MEq->Close();
 260     if( lCalibFile_V0AEq ) lCalibFile_V0AEq->Close();
 261     if( lCalibFile_V0CEq ) lCalibFile_V0CEq->Close();
 262
 263     fRunNumber = lLoadThisCalibration; //Loaded!
 264     AliInfo(Form("Finished loading calibration for run %i",lLoadThisCalibration));
 265     return kTRUE;
 266 }