1 //**************************************************************************
2 //* This file is property of and copyright by the ALICE HLT Project *
3 //* ALICE Experiment at CERN, All rights reserved. *
5 //* Primary Authors: leonidas.xaplanteris@gmail.com *
6 //* for The ALICE HLT Project. *
8 //* Permission to use, copy, modify and distribute this software and its *
9 //* documentation strictly for non-commercial purposes is hereby granted *
10 //* without fee, provided that the above copyright notice appears in all *
11 //* copies and that both the copyright notice and this permission notice *
12 //* appear in the supporting documentation. The authors make no claims *
13 //* about the suitability of this software for any purpose. It is *
14 //* provided "as is" without express or implied warranty. *
15 //**************************************************************************
17 #include "AliHLTTriggerFastJet.h"
18 #include "AliESDEvent.h"
19 #include "AliESDtrackCuts.h"
20 #include "AliESDCaloCluster.h"
21 #include "AliHLTCaloClusterReader.h"
22 #include "AliHLTCaloClusterDataStruct.h"
23 #include "AliESDVZERO.h"
24 #include "AliHLTTriggerDecision.h"
25 #include "AliHLTDomainEntry.h"
26 #include "TLorentzVector.h"
27 #include "TRefArray.h"
30 #include "TObjArray.h"
32 #include "AliHLTScalars.h"
34 #include "AliFJWrapper.h"
37 ClassImp(AliHLTTriggerFastJet)
39 AliHLTTriggerFastJet::AliHLTTriggerFastJet() :
43 fFastJetWrapper(NULL),
46 fOCDBEntry("HLT/ConfigHLT/EmcalJetTrigger")
48 fFastJetWrapper = new AliFJWrapper("FastJet","FastJet");
50 // check if the fMakeStats flag is set to use AliHLTScalars
51 // to Add a new scalar:
52 // scalars.Add(const char *name, const char *description, Double_t value)
54 // if ( fMakeStats ) AliHLTScalars scalars;
57 //_____________________________________________________________
59 AliHLTTriggerFastJet::~AliHLTTriggerFastJet() {
62 delete fFastJetWrapper;
63 fFastJetWrapper = NULL;
70 //_____________________________________________________________
72 int AliHLTTriggerFastJet::DoInit(int argc, const char** argv) {
74 int iResult = ConfigureFromCDBTObjString(fOCDBEntry);
76 if ( iResult>=0 && argc>0 ) {
77 iResult = ConfigureFromArgumentString(argc, argv);
78 HLTImportant("Trigger threshold set from argument string : %.02f GeV:", fEThreshold );
79 } else if ( iResult>=0 ) {
80 HLTImportant("Trigger threshold set from OCDB database entry : %.02f Gev:", fEThreshold );
84 //______________________________________________________________
86 int AliHLTTriggerFastJet::DoDeInit() {
91 //______________________________________________________________
95 AliHLTComponent* AliHLTTriggerFastJet::Spawn() {
96 // see header file for class documentation
97 return new AliHLTTriggerFastJet;
101 const char* AliHLTTriggerFastJet::GetTriggerName() const {
102 // see header file for class documentation
103 return "EmcalJetTrigger";
106 Int_t AliHLTTriggerFastJet::DoTrigger() {
109 AliHLTScalars scalars;
111 if ( GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR ) )
114 const TObject *obj = GetFirstInputObject( kAliHLTAllDataTypes, "AliESDEvent" );
115 AliESDEvent *esd = dynamic_cast<AliESDEvent*>(const_cast<TObject*>(obj));
117 // check for MatchedTrack to avoid double counting
122 esd->GetStdContent();
124 Double_t vertex[3] = {0,0,0};
125 vertex[0] = esd->GetVertex()->GetX();
126 vertex[1] = esd->GetVertex()->GetY();
127 vertex[2] = esd->GetVertex()->GetZ();
128 //cout << "Vertex " << vertex[0] << vertex[1] << vertex [2] << endl;
129 //TLorentzVector gamma;
132 // uncomment for usage (to avoid warning)
135 AliESDVZERO *v0 = esd->GetVZEROData();
137 Float_t MtotVOA = v0->GetMTotV0A();
138 Float_t MtotVOC = v0->GetMTotV0A();
140 Int_t nTracks = esd->GetNumberOfTracks();
142 TObjArray *tracks = new TObjArray(nTracks);
146 EsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
147 EsdTrackCuts->SetMinNClustersTPC(70);
149 // add the tracks to FastJet Wrapper
151 for ( Int_t i = 0; i < nTracks; i++ ) {
153 const AliESDVertex *vtxTPC = esd->GetPrimaryVertexTPC();
154 if (!vtxTPC) continue;
155 //cout << "Got primary vertex " << endl;
157 AliESDtrack *esdtrack = esd->GetTrack(i);
158 if (!esdtrack) continue;
159 tracks->Add(esdtrack);
160 //cout << "Got esd track " << endl;
162 AliESDtrack *tpctrack = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(esd),esdtrack->GetID());
163 if (!tpctrack) continue;
164 //cout << "Got tpc track " << endl;
166 AliExternalTrackParam exParam;
167 Bool_t relate = tpctrack->RelateToVertexTPC(vtxTPC,esd->GetMagneticField(),kVeryBig,&exParam);
173 tpctrack->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
174 //cout << "Set param to tpctrack " << endl;
176 fFastJetWrapper->AddInputVector(tpctrack->Px(),tpctrack->Py(),tpctrack->Pz(),tpctrack->P());
177 if ( fMakeStats ) scalars.Add("TracksPt","TPC tracks pT", tpctrack->Pt());
178 //cout << "Added track with P: " << tpctrack->P() << endl;
182 // // add the Calo Clusters to FastJet Wrapper
183 // TRefArray *caloClustersArr = new TRefArray();
184 // esd->GetEMCALClusters(caloClustersArr);
185 // const Int_t nEMCALClusters = caloClustersArr->GetEntries();
187 // for ( UInt_t ic = 0; ic < nEMCALClusters; ic++ )
189 // AliESDCaloCluster *cluster = dynamic_cast<AliESDCaloCluster*>( caloClustersArr->At(ic) );
190 // Double_t clusterEn = cluster->E();
191 // Double_t subE = 0;
193 // TArrayI *matched = cluster->GetTracksMatched();
194 // for ( UInt_t jk = 0; jk < matched->GetSize(); jk++ )
196 // if ( matched->At(jk) > -1 && matched->At(jk) < esd ->GetNumberOfTracks() )
198 // AliESDtrack *track = dynamic_cast<AliESDtrack*>( tracks->At( (matched->At(jk)) ) );
199 // if ( !track ) continue;
200 // subE += track->P();
203 // clusterEn -= subE;
204 // if ( clusterEn < 0)
206 // cluster->SetE( clusterEn );
207 // if ( cluster->E() == 0 ) continue;
209 // cluster->GetMomentum( gamma, vertex );
210 // fFastJetWrapper->AddInputVector( gamma.Px(), gamma.Py(), gamma.Pz(), clusterEn );
211 // if ( fMakeStats ) {
212 // scalars.Add("ClusterEn","Clusters Energy",clusterEn);
213 // scalars.Add("ClusterEta","Clusters Eta", gamma.Eta());
214 // scalars.Add("ClusterPhi","Clusters Phi", gamma.Phi());
219 for ( Int_t j = 0; j< esd->GetNumberOfCaloClusters(); j++) {
220 //AliESDCaloCluster *cluster = static_cast<AliESDCaloCluster*>(esd->GetCaloCluster(j));
221 AliVCluster *cluster = dynamic_cast<AliVCluster*>(esd->GetCaloCluster(j));
222 if ( cluster->IsEMCAL() ) {
223 Int_t trackindex = cluster->GetTrackMatchedIndex();
224 if (trackindex<0) continue;
225 AliESDtrack *track = esd->GetTrack(trackindex);
226 if (!track) continue;
227 fFastJetWrapper->AddInputVector(track->Px(),track->Py(),track->Pz(),cluster->E());
228 //cout << "Added cluster with E: " << cluster->E() << endl;
230 scalars.Add("ClusterEn" , "Clusters Energy", cluster->E());
231 scalars.Add("ClusterEta", "Clusters Eta" , track->Eta());
232 scalars.Add("ClusterPhi", "Clusters Phi" , track->Phi());
238 fFastJetWrapper->SetupStrategyfromOpt("Best");
239 fFastJetWrapper->SetupAlgorithmfromOpt("antikt");
240 fFastJetWrapper->SetupSchemefromOpt("BIpt");
241 fFastJetWrapper->SetupAreaTypefromOpt("active");
242 fFastJetWrapper->SetNRepeats(1);
243 fFastJetWrapper->SetGhostArea(0.01);
244 fFastJetWrapper->SetMaxRap(0.9);
245 fFastJetWrapper->SetR(0.4);
250 fFastJetWrapper->Run();
251 fFastJetWrapper->GetMedianAndSigma(median,sigma);
253 // checks if trigger criteria are fulfilled
254 if ( TriggerOnJet(fFastJetWrapper->GetSubtractedJetsPts(median)) ) {
255 if ( fMakeStats ) iResult = PushBack(&scalars, kAliHLTDataTypeEventStatistics|kAliHLTDataOriginHLT);
259 // if we got to this point it meens the trigger criteria were not fulfilled
261 description.Form(" No jets with energy > %.02f GeV found! ", fDetector.Data(), fEThreshold);
262 SetDescription(description.Data());
263 TriggerEvent(kFALSE);
265 if ( fMakeStats ) iResult = PushBack(&scalars, kAliHLTDataTypeEventStatistics|kAliHLTDataOriginHLT);
269 //______________________________________________________________
272 Bool_t AliHLTTriggerFastJet::TriggerOnJet(T Jet) {
274 AliHLTScalars scalars;
276 for ( unsigned int ij = 0; ij<Jet.size(); ij++ ) {
277 if ( Jet[ij] > fEThreshold ) {
279 description.Form(" Event contains at least one %s jet with energy > %.02f GeV! ", fDetector.Data(), fEThreshold);
280 SetDescription(description.Data());
282 // Enable the detectors for readout
283 GetReadoutList().Enable(AliHLTReadoutList::kEMCAL |
284 AliHLTReadoutList::kTPC);
286 // Add the available HLT info for readout
287 GetTriggerDomain().Add(kAliHLTAnyDataTypeID, fDetector.Data());
289 // Set trigger desicion
292 // if fMakeStats true fill stats
293 if ( fMakeStats ) scalars.Add("JetsPt","Jets pT", Jet[ij]);
302 //______________________________________________________________
304 int AliHLTTriggerFastJet::Reconfigure(const char* cdbEntry, const char* /*chainId*/) {
306 const char* entry = cdbEntry;
307 if ( !entry || entry[0]==0 ) entry=fOCDBEntry;
309 return ConfigureFromCDBTObjString(entry);
313 //______________________________________________________________
315 int AliHLTTriggerFastJet::ScanConfigurationArgument(int argc, const char** argv) {
317 if ( argc<=0 ) return 0;
319 TString argument = argv[i];
321 if ( argument.CompareTo("-energy") == 0 ) {
322 if (++i>=argc) return -EPROTO;
324 fEThreshold=argument.Atof();
328 if ( argument.CompareTo("-makestats") == 0) {
336 //______________________________________________________________
338 void AliHLTTriggerFastJet::GetOutputDataTypes(AliHLTComponentDataTypeList &list) const {
339 // return the output data types generated
341 list.push_back(kAliHLTDataTypeTriggerDecision);
342 list.push_back(kAliHLTDataTypeEventStatistics|kAliHLTDataOriginHLT);
344 //_______________________________________________________________
346 void AliHLTTriggerFastJet::GetOutputDataSize(unsigned long &constBase, double &inputMultiplier) {
348 constBase = sizeof(AliHLTTriggerDecision) + sizeof(AliHLTDomainEntry)*14;
352 //______________________________________________________________
354 void AliHLTTriggerFastJet::GetOCDBObjectDescription(TMap* const targetMap) {
356 if ( !targetMap ) return;
357 targetMap->Add(new TObjString(fOCDBEntry),
358 new TObjString(Form("%s threshold trigger OCDB object", fDetector.Data()) )
362 //_______________________________________________________________