1 #include "AliAnalysisEtSelectorPhos.h"
2 #include "AliAnalysisEtCuts.h"
3 #include "AliESDCaloCluster.h"
6 #include "AliPHOSGeometry.h"
12 ClassImp(AliAnalysisEtSelectorPhos)
13 AliAnalysisEtSelectorPhos::AliAnalysisEtSelectorPhos(AliAnalysisEtCuts* cuts): AliAnalysisEtSelector(cuts)
22 AliAnalysisEtSelectorPhos::~AliAnalysisEtSelectorPhos()
27 TRefArray* AliAnalysisEtSelectorPhos::GetClusters()
32 Int_t AliAnalysisEtSelectorPhos::Init(Int_t runNumber)
34 AliAnalysisEtSelector::Init(runNumber);
36 int res = LoadGeometry();
41 Bool_t AliAnalysisEtSelectorPhos::CutMinEnergy(const AliESDCaloCluster& cluster) const
43 return cluster.E() > fCuts->GetReconstructedPhosClusterEnergyCut();
46 Bool_t AliAnalysisEtSelectorPhos::CutDistanceToBadChannel(const AliESDCaloCluster& cluster) const
49 cluster.GetPosition(gPos);
52 fGeoUtils->GlobalPos2RelId(glVec, relId);
55 fGeoUtils->Global2Local(locVec, glVec, relId[0]);
56 // std::cout << fGeoUtils << std::endl;
57 //std::cout << relId[0] << " " << cluster.IsPHOS() << std::endl;
58 //std::cout << locVec[0] << " " << " " << locVec[1] << " " << locVec[2] << std::endl;
59 for (Int_t x = 0; x < fBadMapM2->GetNbinsX(); x++)
61 for (Int_t z = 0; z < fBadMapM2->GetNbinsY(); z++)
65 if (fBadMapM2->GetBinContent(x+1, z+1) != 0)
75 fGeoUtils->RelPosInModule(tmpRel, tmpX, tmpZ);
77 Float_t distance = TMath::Sqrt((tmpX-locVec[0])*(tmpX-locVec[0]) + (tmpZ - locVec[2])*(tmpZ-locVec[2]));
78 //Float_t distance = TMath::Sqrt((x-relId[3])*(x-relId[3]) + (z - relId[2])*(z-relId[2]));
80 if (distance < fCuts->GetPhosBadDistanceCut())
82 // std::cout << "Module 2, position: " << locVec[0] << ", " << locVec[2] << ", distance to bad channel: " << distance << ", number of cells: " << cluster.GetNCells() << std::endl;
89 if (fBadMapM3->GetBinContent(x+1, z+1) != 0)
99 fGeoUtils->RelPosInModule(tmpRel, tmpX, tmpZ);
101 Float_t distance = TMath::Sqrt((tmpX-locVec[0])*(tmpX-locVec[0]) + (tmpZ - locVec[2])*(tmpZ-locVec[2]));
103 // Float_t distance = TMath::Sqrt((x-locVec[0])*(x-locVec[0]) + (z - locVec[2])*(z-locVec[2]));
104 //Float_t distance = TMath::Sqrt((x-relId[3])*(x-relId[3]) + (z - relId[2])*(z-relId[2]));
105 if (distance < fCuts->GetPhosBadDistanceCut())
107 // std::cout << "Module 3, position: " << locVec[0] << ", " << locVec[2] << ", distance to bad channel: " << distance << ", number of cells: " << cluster.GetNCells() << std::endl;
114 if (fBadMapM4->GetBinContent(x+1, z+1) != 0)
124 fGeoUtils->RelPosInModule(tmpRel, tmpX, tmpZ);
126 Float_t distance = TMath::Sqrt((tmpX-locVec[0])*(tmpX-locVec[0]) + (tmpZ - locVec[2])*(tmpZ-locVec[2]));
128 // Float_t distance = TMath::Sqrt((x-locVec[0])*(x-locVec[0]) + (z - locVec[2])*(z-locVec[2]));
129 //Float_t distance = TMath::Sqrt((x-relId[3])*(x-relId[3]) + (z - relId[2])*(z-relId[2]));
130 if (distance < fCuts->GetPhosBadDistanceCut())
132 // std::cout << "Module 4, position: " << locVec[0] << ", " << locVec[2] << ", distance to bad channel: " << distance << ", number of cells: " << cluster.GetNCells() << std::endl;
145 Bool_t AliAnalysisEtSelectorPhos::CutTrackMatching(const AliESDCaloCluster& cluster, Double_t &r) const
149 // cluster->GetTrackDx(), cluster->GetTrackDz(), event->GetTrack(trackMatchedIndex)->Pt(), event->GetTrack(trackMatchedIndex)->Charge(), ev
151 Int_t nTracksMatched = cluster.GetNTracksMatched();
152 if(nTracksMatched == 0) return kFALSE;
154 Int_t trackMatchedIndex = cluster.GetTrackMatchedIndex();
155 if(trackMatchedIndex < 0) return kTRUE;
157 AliVParticle *track = fEvent->GetTrack(trackMatchedIndex);
158 Double_t dx = cluster.GetTrackDx();
159 Double_t dz = cluster.GetTrackDz();
160 Double_t pt = track->Pt();
161 Int_t charge = track->Charge();
165 Double_t sx=TMath::Min(5.4,2.59719e+02*TMath::Exp(-pt/1.02053e-01)+
166 6.58365e-01*5.91917e-01*5.91917e-01/((pt-9.61306e-01)*(pt-9.61306e-01)+5.91917e-01*5.91917e-01)+1.59219);
167 Double_t sz=TMath::Min(2.75,4.90341e+02*1.91456e-02*1.91456e-02/(pt*pt+1.91456e-02*1.91456e-02)+1.60) ;
169 Double_t mf = fEvent->GetMagneticField(); //Positive for ++ and negative for --
171 if(mf<0.){ //field --
174 meanX=TMath::Min(7.3, 3.89994*1.20679*1.20679/(pt*pt+1.20679*1.20679)+0.249029+2.49088e+07*TMath::Exp(-pt*3.33650e+01)) ;
176 meanX=-TMath::Min(7.7,3.86040*0.912499*0.912499/(pt*pt+0.912499*0.912499)+1.23114+4.48277e+05*TMath::Exp(-pt*2.57070e+01)) ;
181 meanX=-TMath::Min(8.0,3.86040*1.31357*1.31357/(pt*pt+1.31357*1.31357)+0.880579+7.56199e+06*TMath::Exp(-pt*3.08451e+01)) ;
183 meanX= TMath::Min(6.85, 3.89994*1.16240*1.16240/(pt*pt+1.16240*1.16240)-0.120787+2.20275e+05*TMath::Exp(-pt*2.40913e+01)) ;
186 Double_t rz=(dz-meanZ)/sz ;
187 Double_t rx=(dx-meanX)/sx ;
188 r = TMath::Sqrt(rx*rx+rz*rz);
189 if(r < fCuts->GetPhosTrackRCut()) return kFALSE;
195 int AliAnalysisEtSelectorPhos::LoadGeometry()
198 fGeoUtils = AliPHOSGeometry::GetInstance("IHEP");
199 // ifstream f("badchannels.txt", ios::in);
203 int AliAnalysisEtSelectorPhos::LoadBadMaps()
205 TFile *f = TFile::Open("badchannels.root", "READ");
209 std::cout << "Could not open badchannels.root" << std::endl;
213 fBadMapM2 = (TH2I*)f->Get("bad_channels_m2");
216 std::cout << "Could not find bad_channels_m2 in badchannels.root" << std::endl;
218 fBadMapM3 = (TH2I*)f->Get("bad_channels_m3");
221 std::cout << "Could not find bad_channels_m3 in badchannels.root" << std::endl;
224 fBadMapM4 = (TH2I*)f->Get("bad_channels_m4");
227 std::cout << "Could not find bad_channels_m4 in badchannels.root" << std::endl;