1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 //-------------------------------------------------------------------------
18 // AOD reader for jet analysis
19 // This is the reader which must be used if the jet analysis task
20 // is executed after the ESD filter task, in order to read its output
22 // Author: Davide Perrino <davide.perrino@cern.ch>
23 //-------------------------------------------------------------------------
26 #include <Riostream.h>
28 #include <TLorentzVector.h>
33 #include <TGeoManager.h>
35 #include "AliJetAODReader.h"
36 #include "AliJetAODReaderHeader.h"
37 #include "AliAODEvent.h"
38 #include "AliAODTrack.h"
39 #include "AliAODMCParticle.h"
40 #include "AliJetDummyGeo.h"
41 #include "AliJetAODFillUnitArrayTracks.h"
42 #include "AliJetAODFillUnitArrayEMCalDigits.h"
43 #include "AliJetHadronCorrectionv1.h"
44 #include "AliJetUnitArray.h"
46 ClassImp(AliJetAODReader)
48 AliJetAODReader::AliJetAODReader():
64 fApplyElectronCorrection(kFALSE),
65 fApplyMIPCorrection(kTRUE),
66 fApplyFractionHadronicCorrection(kFALSE),
67 fFractionHadronicCorrection(0.3),
81 //____________________________________________________________________________
83 AliJetAODReader::~AliJetAODReader()
101 //____________________________________________________________________________
103 void AliJetAODReader::OpenInputFiles()
105 // Open the necessary input files
106 // chain for the AODs
107 fChain = new TChain("aodTree");
109 // get directory and pattern name from the header
110 const char* dirName=fReaderHeader->GetDirectory();
111 const char* pattern=fReaderHeader->GetPattern();
113 // // Add files matching patters to the chain
115 void *dir = gSystem->OpenDirectory(dirName);
116 const char *name = 0x0;
117 int naod = ((AliJetAODReaderHeader*) fReaderHeader)->GetNaod();
119 while ((name = gSystem->GetDirEntry(dir))){
120 if (a>=naod) continue;
122 if (strstr(name,pattern)){
124 sprintf(path,"%s/%s/aod.root",dirName,name);
125 fChain->AddFile(path);
130 gSystem->FreeDirectory(dir);
133 fChain->SetBranchAddress("AOD",&fAOD);
135 int nMax = fChain->GetEntries();
137 printf("\n AliJetAODReader: Total number of events in chain= %d \n",nMax);
139 // set number of events in header
140 if (fReaderHeader->GetLastEvent() == -1)
141 fReaderHeader->SetLastEvent(nMax);
143 Int_t nUsr = fReaderHeader->GetLastEvent();
144 fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
148 //____________________________________________________________________________
150 void AliJetAODReader::ConnectTree(TTree* tree, TObject* /*data*/) {
152 // For AOD reader it's needed only to set the number of events
153 fChain = (TChain*) tree;
155 Int_t nMax = fChain->GetEntries();
156 printf("\n AliJetAODReader: Total number of events in chain= %5d \n", nMax);
157 // set number of events in header
158 if (fReaderHeader->GetLastEvent() == -1)
159 fReaderHeader->SetLastEvent(nMax);
161 Int_t nUsr = fReaderHeader->GetLastEvent();
162 fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
167 Bool_t AliJetAODReader::AcceptAODMCParticle(AliAODMCParticle *mcP,Short_t flag){
170 // filter for charge and for charged and neutral, no detector
172 // physical priamries already selected
174 Int_t pdg = TMath::Abs(mcP->GetPdgCode());
176 // exclude neutrinos anyway
177 if((pdg == 12 || pdg == 14 || pdg == 16)) return kFALSE;
179 if(flag==AliJetAODReaderHeader::kAllMC)return kTRUE;
180 if(flag==AliJetAODReaderHeader::kChargedMC){
181 if(mcP->Charge()==0)return kFALSE;
190 Bool_t AliJetAODReader::FillMomentumArrayMC(){
193 // This routine fetches the MC particles from the AOD
194 // Depending on the flag all particles except neurinos are use
195 // or only charged particles
198 TClonesArray *mcArray = (TClonesArray*)fAOD->FindListObject(AliAODMCParticle::StdBranchName());
200 Printf("%s:%d No MC particle branch found",(char*)__FILE__,__LINE__);
204 Int_t nMC = mcArray->GetEntriesFast();
206 // get number of tracks in event (for the loop)
207 if(fDebug>0)printf("AOD MC tracks: %5d \t", nMC);
209 // temporary storage of signal and pt cut flag
210 Int_t* sflag = new Int_t[nMC];
211 Int_t* cflag = new Int_t[nMC];
213 // get cuts set by user
214 Float_t ptMin = fReaderHeader->GetPtCut();
215 Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
216 Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
221 Short_t readerFlag = ((AliJetAODReaderHeader*)fReaderHeader)->GetReadAODMC();
224 for (Int_t it = 0; it < nMC; it++) {
225 AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
226 if(!track->IsPhysicalPrimary())continue;
228 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
230 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
231 if(!AcceptAODMCParticle(track,readerFlag))continue;
237 if (mcTrack == 0) new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
238 new ((*fMomentumArray)[mcTrack]) TLorentzVector(p3,p3.Mag());
240 cflag[mcTrack] = ( pt > ptMin ) ? 1: 0;
244 if(fDebug>0)printf("Used MC tracks: %5d \n", mcTrack);
245 // set the signal flags
246 fSignalFlag.Set(mcTrack,sflag);
247 fCutFlag.Set(mcTrack,cflag);
255 //____________________________________________________________________________
257 Bool_t AliJetAODReader::FillMomentumArray()
259 // Clear momentum array
262 fDebug = fReaderHeader->GetDebug();
264 if(((AliJetAODReaderHeader*)fReaderHeader)->GetReadAODMC()){
265 return FillMomentumArrayMC();
274 // get number of tracks in event (for the loop)
275 Int_t nt = fAOD->GetNTracks();
276 if(fDebug>0)printf("AOD tracks: %5d \t", nt);
278 // temporary storage of signal and pt cut flag
279 Int_t* sflag = new Int_t[nt];
280 Int_t* cflag = new Int_t[nt];
282 // get cuts set by user
283 Float_t ptMin = fReaderHeader->GetPtCut();
284 Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
285 Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
286 UInt_t filterMask = ((AliJetAODReaderHeader*)fReaderHeader)->GetTestFilterMask();
293 for (Int_t it = 0; it < nt; it++) {
294 AliAODTrack *track = fAOD->GetTrack(it);
295 UInt_t status = track->GetStatus();
297 Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
298 p3.SetXYZ(mom[0],mom[1],mom[2]);
301 if (status == 0) continue;
302 if((filterMask>0)&&!(track->TestFilterBit(filterMask)))continue;
303 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
305 if (aodTrack == 0) new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
307 new ((*fMomentumArray)[aodTrack]) TLorentzVector(p3,p3.Mag());
308 sflag[aodTrack] = (TMath::Abs(track->GetLabel()) < 10000) ? 1 : 0;
309 cflag[aodTrack] = ( pt > ptMin ) ? 1: 0;
313 if(fDebug>0)printf("Used AOD tracks: %5d \n", aodTrack);
314 // set the signal flags
315 fSignalFlag.Set(aodTrack,sflag);
316 fCutFlag.Set(aodTrack,cflag);
324 //__________________________________________________________
325 void AliJetAODReader::SetApplyMIPCorrection(Bool_t val)
328 // Set flag to apply MIP correction fApplyMIPCorrection
329 // - exclusive with fApplyFractionHadronicCorrection
332 fApplyMIPCorrection = val;
333 if(fApplyMIPCorrection == kTRUE)
335 SetApplyFractionHadronicCorrection(kFALSE);
336 printf("Enabling MIP Correction \n");
340 printf("Disabling MIP Correction \n");
344 //__________________________________________________________
345 void AliJetAODReader::SetApplyFractionHadronicCorrection(Bool_t val)
348 // Set flag to apply EMC hadronic correction fApplyFractionHadronicCorrection
349 // - exclusive with fApplyMIPCorrection
352 fApplyFractionHadronicCorrection = val;
353 if(fApplyFractionHadronicCorrection == kTRUE)
355 SetApplyMIPCorrection(kFALSE);
356 printf("Enabling Fraction Hadronic Correction \n");
360 printf("Disabling Fraction Hadronic Correction \n");
364 //__________________________________________________________
365 void AliJetAODReader::SetFractionHadronicCorrection(Double_t val)
368 // Set value to fFractionHadronicCorrection (default is 0.3)
369 // apply EMC hadronic correction fApplyFractionHadronicCorrection
370 // - exclusive with fApplyMIPCorrection
373 fFractionHadronicCorrection = val;
374 if(fFractionHadronicCorrection > 0.0 && fFractionHadronicCorrection <= 1.0)
376 SetApplyFractionHadronicCorrection(kTRUE);
377 printf("Fraction Hadronic Correction %1.3f \n",fFractionHadronicCorrection);
381 SetApplyFractionHadronicCorrection(kFALSE);
385 //____________________________________________________________________________
386 void AliJetAODReader::CreateTasks(TChain* tree)
389 // For reader task initialization
392 fDebug = fReaderHeader->GetDebug();
393 fDZ = fReaderHeader->GetDZ();
396 // Init EMCAL geometry and create UnitArray object
398 // cout << "In create task" << endl;
402 fFillUnitArray = new TTask("fFillUnitArray","Fill unit array jet finder");
403 fFillUAFromTracks = new AliJetAODFillUnitArrayTracks();
404 fFillUAFromTracks->SetReaderHeader(fReaderHeader);
405 fFillUAFromTracks->SetGeom(fGeom);
406 fFillUAFromTracks->SetTPCGrid(fTpcGrid);
407 fFillUAFromTracks->SetEMCalGrid(fEmcalGrid);
411 fFillUAFromTracks->SetGrid0(fGrid0);
412 fFillUAFromTracks->SetGrid1(fGrid1);
413 fFillUAFromTracks->SetGrid2(fGrid2);
414 fFillUAFromTracks->SetGrid3(fGrid3);
415 fFillUAFromTracks->SetGrid4(fGrid4);
417 fFillUAFromTracks->SetApplyMIPCorrection(fApplyMIPCorrection);
418 fFillUAFromTracks->SetHadCorrector(fHadCorr);
419 fFillUAFromEMCalDigits = new AliJetAODFillUnitArrayEMCalDigits();
420 fFillUAFromEMCalDigits->SetReaderHeader(fReaderHeader);
421 fFillUAFromEMCalDigits->SetGeom(fGeom);
422 fFillUAFromEMCalDigits->SetTPCGrid(fTpcGrid);
423 fFillUAFromEMCalDigits->SetEMCalGrid(fEmcalGrid);
424 fFillUAFromEMCalDigits->SetApplyFractionHadronicCorrection(fApplyFractionHadronicCorrection);
425 fFillUAFromEMCalDigits->SetFractionHadronicCorrection(fFractionHadronicCorrection);
426 fFillUAFromEMCalDigits->SetApplyElectronCorrection(fApplyElectronCorrection);
427 // Add the task to global task
428 fFillUnitArray->Add(fFillUAFromTracks);
429 fFillUnitArray->Add(fFillUAFromEMCalDigits);
430 fFillUAFromTracks->SetActive(kFALSE);
431 fFillUAFromEMCalDigits->SetActive(kFALSE);
433 cout << "Tasks instantiated at that stage ! " << endl;
434 cout << "You can loop over events now ! " << endl;
438 //____________________________________________________________________________
439 Bool_t AliJetAODReader::ExecTasks(Bool_t procid, TRefArray* refArray)
443 // Fill the reader part
447 fRefArray = refArray;
448 //(not used ?) Int_t nEntRef = fRefArray->GetEntries();
449 //(not used ?) Int_t nEntUnit = fUnitArray->GetEntries();
451 // clear momentum array
455 fDebug = fReaderHeader->GetDebug();
456 fOpt = fReaderHeader->GetDetector();
462 // TPC only or Digits+TPC or Clusters+TPC
463 if(fOpt%2==!0 && fOpt!=0){
464 fFillUAFromTracks->SetAOD(fAOD);
465 fFillUAFromTracks->SetActive(kTRUE);
466 fFillUAFromTracks->SetUnitArray(fUnitArray);
467 fFillUAFromTracks->SetRefArray(fRefArray);
468 fFillUAFromTracks->SetReferences(fRef);
469 fFillUAFromTracks->SetSignalFlag(fSignalFlag);
470 fFillUAFromTracks->SetCutFlag(fCutFlag);
471 fFillUAFromTracks->SetProcId(fProcId);
472 // fFillUAFromTracks->ExecuteTask("tpc"); // => Temporarily changed
473 fFillUAFromTracks->Exec("tpc");
475 fNumCandidate = fFillUAFromTracks->GetMult();
476 fNumCandidateCut = fFillUAFromTracks->GetMultCut();
478 fSignalFlag = fFillUAFromTracks->GetSignalFlag();
479 fCutFlag = fFillUAFromTracks->GetCutFlag();
482 // Digits only or Digits+TPC
483 if(fOpt>=2 && fOpt<=3){
484 fFillUAFromEMCalDigits->SetAOD(fAOD);
485 fFillUAFromEMCalDigits->SetActive(kTRUE);
486 fFillUAFromEMCalDigits->SetUnitArray(fUnitArray);
487 fFillUAFromEMCalDigits->SetRefArray(fRefArray);
488 fFillUAFromEMCalDigits->SetProcId(fFillUAFromTracks->GetProcId());
489 fFillUAFromEMCalDigits->SetInitMult(fFillUAFromTracks->GetMult());
490 fFillUAFromEMCalDigits->SetInitMultCut(fFillUAFromTracks->GetMultCut());
491 fFillUAFromEMCalDigits->Exec("digits"); // => Temporarily added
492 fNumCandidate = fFillUAFromEMCalDigits->GetMult();
493 fNumCandidateCut = fFillUAFromEMCalDigits->GetMultCut();
496 // fFillUnitArray->ExecuteTask(); // => Temporarily commented
501 //____________________________________________________________________________
502 Bool_t AliJetAODReader::SetEMCALGeometry()
505 // Set the EMCal Geometry
508 if (!fTree->GetFile())
511 TString geomFile(fTree->GetFile()->GetName());
512 geomFile.ReplaceAll("AliESDs", "geometry");
514 // temporary workaround for PROOF bug #18505
515 geomFile.ReplaceAll("#geometry.root#geometry.root", "#geometry.root");
516 if(fDebug>1) printf("Current geometry file %s \n", geomFile.Data());
518 // Define EMCAL geometry to be able to read ESDs
519 fGeom = AliJetDummyGeo::GetInstance();
521 fGeom = AliJetDummyGeo::GetInstance("EMCAL_COMPLETE","EMCAL");
523 // To be setted to run some AliEMCALGeometry functions
524 TGeoManager::Import(geomFile);
525 fGeom->GetTransformationForSM();
526 printf("\n EMCal Geometry set ! \n");
532 //____________________________________________________________________________
533 void AliJetAODReader::InitParameters()
535 // Initialise parameters
536 fOpt = fReaderHeader->GetDetector();
537 fNumUnits = fGeom->GetNCells(); // Number of cells in EMCAL
538 if(fDebug>1) printf("\n EMCal parameters initiated ! \n");
541 //____________________________________________________________________________
542 void AliJetAODReader::InitUnitArray()
544 //Initialises unit arrays
545 Int_t nElements = fTpcGrid->GetNEntries();
546 Float_t eta = 0., phi = 0., deltaEta = 0., deltaPhi = 0.;
547 if(fArrayInitialised) fUnitArray->Delete();
549 if(fTpcGrid->GetGridType()==0)
550 { // Fill the following quantities :
551 // Good track ID, (Eta,Phi) position ID, eta, phi, energy, px, py, pz, deltaEta, deltaPhi,
552 // detector flag, in/out jet, pt cut, mass, cluster ID)
553 for(Int_t nBin = 1; nBin < nElements+1; nBin++)
555 // fTpcGrid->GetEtaPhiFromIndex2(nBin,eta,phi);
556 fTpcGrid->GetEtaPhiFromIndex2(nBin,phi,eta);
557 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
558 deltaEta = fTpcGrid->GetDeta();
559 deltaPhi = fTpcGrid->GetDphi();
560 new ((*fUnitArray)[nBin-1]) AliJetUnitArray(nBin-1,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
564 if(fTpcGrid->GetGridType()==1)
567 Int_t n0 = 0, n1 = 0, n2 = 0, n3 = 0, n4 = 0;
571 // Define a grid of cell for the gaps between SM
572 Double_t phimin0 = 0., phimin1 = 0., phimin2 = 0., phimin3 = 0., phimin4 = 0.;
573 Double_t phimax0 = 0., phimax1 = 0., phimax2 = 0., phimax3 = 0., phimax4 = 0.;
574 fGeom->GetPhiBoundariesOfSMGap(0,phimin0,phimax0);
575 fGrid0 = new AliJetGrid(0,95,phimin0,phimax0,-0.7,0.7); // 0.015 x 0.015
576 fGrid0->SetGridType(0);
577 fGrid0->SetMatrixIndexes();
578 fGrid0->SetIndexIJ();
579 n0 = fGrid0->GetNEntries();
580 fGeom->GetPhiBoundariesOfSMGap(1,phimin1,phimax1);
581 fGrid1 = new AliJetGrid(0,95,phimin1,phimax1,-0.7,0.7); // 0.015 x 0.015
582 fGrid1->SetGridType(0);
583 fGrid1->SetMatrixIndexes();
584 fGrid1->SetIndexIJ();
585 n1 = fGrid1->GetNEntries();
586 fGeom->GetPhiBoundariesOfSMGap(2,phimin2,phimax2);
587 fGrid2 = new AliJetGrid(0,95,phimin2,phimax2,-0.7,0.7); // 0.015 x 0.015
588 fGrid2->SetGridType(0);
589 fGrid2->SetMatrixIndexes();
590 fGrid2->SetIndexIJ();
591 n2 = fGrid2->GetNEntries();
592 fGeom->GetPhiBoundariesOfSMGap(3,phimin3,phimax3);
593 fGrid3 = new AliJetGrid(0,95,phimin3,phimax3,-0.7,0.7); // 0.015 x 0.015
594 fGrid3->SetGridType(0);
595 fGrid3->SetMatrixIndexes();
596 fGrid3->SetIndexIJ();
597 n3 = fGrid3->GetNEntries();
598 fGeom->GetPhiBoundariesOfSMGap(4,phimin4,phimax4);
599 fGrid4 = new AliJetGrid(0,95,phimin4,phimax4,-0.7,0.7); // 0.015 x 0.015
600 fGrid4->SetGridType(0);
601 fGrid4->SetMatrixIndexes();
602 fGrid4->SetIndexIJ();
603 n4 = fGrid4->GetNEntries();
605 nGaps = n0+n1+n2+n3+n4;
609 for(Int_t nBin = 0; nBin < fNumUnits+nElements+nGaps; nBin++)
613 fGeom->EtaPhiFromIndex(nBin, eta, phi); // From EMCal geometry
614 // fEmcalGrid->GetEtaPhiFromIndex2(nBin,phi,eta); // My function from Grid
615 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
616 deltaEta = fEmcalGrid->GetDeta(); // Modify with the exact detector values
617 deltaPhi = fEmcalGrid->GetDphi(); // Modify with the exact detector values
618 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
621 if(nBin>=fNumUnits && nBin<fNumUnits+nElements){
622 fTpcGrid->GetEtaPhiFromIndex2(nBin+1-fNumUnits,phi,eta);
623 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
624 deltaEta = fTpcGrid->GetDeta();
625 deltaPhi = fTpcGrid->GetDphi();
626 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
630 if(nBin>=fNumUnits+nElements && nBin<fNumUnits+nElements+nGaps){
631 if(nBin<fNumUnits+nElements+n0)
634 fGrid0->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements),phi,eta);
635 deltaEta = fGrid0->GetDeta();
636 deltaPhi = fGrid0->GetDphi();
637 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
639 else if(nBin>=fNumUnits+nElements+n0 && nBin<fNumUnits+nElements+n0+n1)
642 fGrid1->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0),phi,eta);
643 deltaEta = fGrid1->GetDeta();
644 deltaPhi = fGrid1->GetDphi();
645 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
647 else if(nBin>=fNumUnits+nElements+n0+n1 && nBin<fNumUnits+nElements+n0+n1+n2)
650 fGrid2->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1),phi,eta);
651 deltaEta = fGrid2->GetDeta();
652 deltaPhi = fGrid2->GetDphi();
653 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
655 else if(nBin>=fNumUnits+nElements+n0+n1+n2 && nBin<fNumUnits+nElements+n0+n1+n2+n3)
658 fGrid3->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1+n2),phi,eta);
659 deltaEta = fGrid3->GetDeta();
660 deltaPhi = fGrid3->GetDphi();
661 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
663 else if(nBin>=fNumUnits+nElements+n0+n1+n2+n3 && nBin<fNumUnits+nElements+nGaps)
666 fGrid4->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1+n2+n3),phi,eta);
667 deltaEta = fGrid4->GetDeta();
668 deltaPhi = fGrid4->GetDphi();
669 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
675 } // end loop on nBin
676 } // end grid type == 1
677 fArrayInitialised = 1;