PWGJE/UserTasks/AliAnalysisTaskJetShape.h

   1 #ifndef ALIANALYSISTASKJETSHAPE_h
   2 #define ALIANALYSISTASKJETSHAPE_h
   3
   4 #include "AliAnalysisTaskSE.h"
   5 #include "AliVEvent.h"
   6 //#include "AliTriggerAnalysis.h"
   7 //#include "THnSparse.h"
   8 #include "TMath.h"
   9 #include "TString.h"
  10 #include "TH1F.h"
  11 #include "TClonesArray.h"
  12
  13 class TH2F;
  14 class TH3F;
  15 class TVector3;
  16
  17 class AliAODExtension;
  18 class TNtuple;
  19 class TTree;
  20 class TList;
  21 class AliAODEvent;
  22 class AliAODJet;
  23 class TArrayD;
  24
  25
  26
  27
  28
  29
  30
  31 class AliAnalysisTaskJetShape : public AliAnalysisTaskSE {
  32
  33  public:
  34
  35
  36 class AliAnalysisTaskJetShapeTool :public TObject {
  37 public:
  38 AliAnalysisTaskJetShapeTool();
  39 //AnalJetSubStrTool(TList *list, TVector3 v);
  40 AliAnalysisTaskJetShapeTool(TClonesArray *list);
  41  ~AliAnalysisTaskJetShapeTool();
  42
  43
  44  void SetVecJ(TVector3 v);
  45  void SetPtMinTr(Double_t a, Double_t b) {fPtMinTr[0] = a; fPtMinTr[1] = b;}
  46
  47  void SetNEntries(Int_t n) {fEntries = n;}
  48
  49  TArrayI GetListJ(Int_t b, Int_t i)  {return fListJ[b][i];}
  50  TArrayI GetListB1(Int_t b, Int_t i) {return fListB1[b][i];}
  51  TArrayI GetListB2(Int_t b, Int_t i) {return fListB2[b][i];}
  52
  53  Int_t GetSizeJ(Int_t b, Int_t i)  {return fListJc[b][i];}
  54  Int_t GetSizeB1(Int_t b, Int_t i) {return fListB1c[b][i];}
  55  Int_t GetSizeB2(Int_t b, Int_t i) {return fListB2c[b][i];}
  56
  57  TVector3 *GetJ(Int_t b,Int_t i, Int_t p)
  58   { return (TVector3*)GetAt(fListJ[b][i].At(p));  }
  59  TVector3 *GetB1(Int_t b,Int_t i, Int_t p)
  60   { return (TVector3*)GetAt(fListB1[b][i].At(p)); }
  61  TVector3 *GetB2(Int_t b,Int_t i, Int_t p)
  62   { return (TVector3*)GetAt(fListB2[b][i].At(p)); }
  63
  64  TVector3 GetVecJ() {return  fvecJ;}
  65  TVector3 GetVecB1() {return fvecB1;}
  66  TVector3 GetVecB2() {return fvecB2;}
  67
  68  Double_t GetR(Int_t b) {return fR[b];}
  69
  70  TVector3 *GetAt(Int_t i)
  71  {
  72    if(i<0 || i>= fEntries) printf(" TVector3 *GetAt(Int_t i) for i= %d\n",i);
  73 return (TVector3*) fList->At(i);
  74 }
  75
  76  Double_t GetLocPhi(TVector3 v, Int_t i)
  77  {
  78    TVector3 p(*GetAt(i));
  79    //   GetLocalMom(v, &vtmp);
  80    //   return CalcdPhi0To2pi(vtmp.Phi());
  81
  82    Double_t phi =  TMath::ATan2(p.Eta() - v.Eta(), CalcdPhi(p.Phi(), v.Phi()) );
  83
  84    return CalcdPhi0To2pi(phi);
  85  }
  86
  87
  88  Double_t GetLocPhiJ(Int_t b, Int_t i, Int_t index)
  89  {
  90    TVector3 v1(*GetAt(fListJ[b][i].At(index)));
  91    Double_t phi =  TMath::ATan2(v1.Eta() - fvecJ.Eta(), CalcdPhi(v1.Phi(), fvecJ.Phi()) );
  92    return CalcdPhi0To2pi(phi);
  93  }
  94
  95
  96  Double_t GetLocPhiB1(Int_t b, Int_t i, Int_t index)
  97  {
  98    TVector3 v0(*GetAt(fListB1[b][i].At(index)));
  99    Double_t cV = (fvecB1(0)*fvecJ(0) + fvecB1(1)*fvecJ(1))/fvecJ.Perp2();
 100    Double_t sV = (fvecB1(1)*fvecJ(0) - fvecB1(0)*fvecJ(1))/fvecJ.Perp2();
 101    TVector3 v1(v0(0)*cV+v0(1)*sV, v0(1)*cV-v0(0)*sV, v0(2));
 102    //   GetLocalMom(fvecJ, &v1);
 103    //   return CalcdPhi0To2pi(v1.Phi());
 104
 105    Double_t phi =  TMath::ATan2(v1.Eta() - fvecJ.Eta(), CalcdPhi(v1.Phi(), fvecJ.Phi()) );
 106    return CalcdPhi0To2pi(phi);
 107
 108  }
 109  Double_t GetLocPhiB2(Int_t b, Int_t i, Int_t index)
 110  {
 111    TVector3 v0(*GetAt(fListB2[b][i].At(index)));
 112    Double_t cV = (fvecB2(0)*fvecJ(0) + fvecB2(1)*fvecJ(1))/fvecJ.Perp2();
 113    Double_t sV = (fvecB2(1)*fvecJ(0) - fvecB2(0)*fvecJ(1))/fvecJ.Perp2();
 114    TVector3 v1(v0(0)*cV+v0(1)*sV, v0(1)*cV-v0(0)*sV, v0(2));
 115    //   GetLocalMom(fvecJ, &v1);
 116    //   return CalcdPhi0To2pi(v1.Phi());
 117    Double_t phi =  TMath::ATan2(v1.Eta() - fvecJ.Eta(), CalcdPhi(v1.Phi(), fvecJ.Phi()) );
 118    return CalcdPhi0To2pi(phi);
 119
 120  }
 121
 122
 123  void Add(TVector3 v)
 124    // { new(fList[fList.GetEntriesFast()]) TVector3(v);}
 125  { new((*fList)[fEntries]) TVector3(v); fEntries++;}
 126
 127  void AddToJ(Int_t b, Int_t i, Int_t index) {
 128    fListJ[b][i].AddAt(index, fListJc[b][i]);
 129    fListJc[b][i]++;
 130  }
 131
 132  void AddToB1(Int_t b, Int_t i, Int_t index) {
 133    fListB1[b][i].AddAt(index, fListB1c[b][i]);
 134    fListB1c[b][i]++;
 135  }
 136
 137  void AddToB2(Int_t b, Int_t i, Int_t index) {
 138    fListB2[b][i].AddAt(index, fListB2c[b][i]);
 139    fListB2c[b][i]++;
 140  }
 141
 142
 143
 144
 145  TVector3 GetPRecJ(Int_t b, Int_t i) {return fPRecJ[b][i];}
 146  TVector3 GetPRecInRJ() {return fPRecInRJ;}//
 147
 148
 149  void Init();
 150  // Bool_t FindCorrelationAxesAnd(TArrayI list, TVector3 vec, Double_t &Phi, Int_t scenario=0);
 151  Bool_t FindCorrelationAxes(TArrayI list, TVector3 vec, Double_t &Phi, Int_t scenario=3);
 152  Bool_t FindCorrelationAxesCorr(TArrayI list, TVector3 vec, Double_t &Phi, Int_t scenario, Double_t R);
 153
 154  // virtual  void  Print(Option_t* /*option = ""*/) const;
 155  void  PRINT() const;
 156  void PRINT(TArrayI a, Int_t n, const char *txt="") const;
 157
 158  void Clean();
 159
 160 static Double_t CalcdPhi0To2pi(Double_t phi1, Double_t phi2=0.)
 161 {
 162   Double_t dphi = CalcdPhi(phi1, phi2);
 163   while(dphi<0) dphi+=TMath::TwoPi();
 164   while(dphi>TMath::TwoPi()) dphi-=TMath::TwoPi();
 165   return dphi;
 166 }
 167
 168 static Double_t CalcdPhi(Double_t phi1, Double_t phi2);
 169
 170 static Double_t CalcdPhiSigned(Double_t phi1, Double_t phi2)
 171 {
 172
 173   Double_t dphi = CalcdPhi(phi1, phi2);
 174   if(dphi <  TMath::Pi()) return dphi;
 175   else return  (dphi-  TMath::Pi());
 176
 177   return -999;
 178 }
 179
 180
 181 private:
 182
 183  TVector3 fvecJ;
 184  TVector3 fvecB1;
 185  TVector3 fvecB2;
 186
 187
 188  Double_t fRmax;
 189
 190  static const Int_t fgkbinR= 3; // n bins
 191  Double_t fR[fgkbinR];
 192  TArrayI fListJ[fgkbinR][2]; //
 193  TArrayI fListB1[fgkbinR][2];//
 194  TArrayI fListB2[fgkbinR][2];//
 195
 196  Int_t fListJc[fgkbinR][2]; //
 197  Int_t fListB1c[fgkbinR][2];//
 198  Int_t fListB2c[fgkbinR][2];//
 199
 200  TVector3 fPRecJ[fgkbinR][2];//
 201  TVector3 fPRecInRJ;//
 202
 203  Double_t fPtMinTr[2];
 204
 205
 206  TClonesArray *fList; //!
 207
 208  Int_t fEntries;
 209
 210  Double_t CalcR(TVector3 v1, TVector3 v2);
 211  void GetLocalMom(TVector3 vecJ, TVector3 *q);
 212
 213
 214   AliAnalysisTaskJetShapeTool(const AliAnalysisTaskJetShapeTool&);            // not implemented
 215   AliAnalysisTaskJetShapeTool& operator=(const AliAnalysisTaskJetShapeTool&); // not implemented
 216
 217  ClassDef(AliAnalysisTaskJetShapeTool, 1)   // tbd
 218
 219 };
 220
 221
 222
 223
 224 class AliAnalysisTaskJetShapeHM :public TObject {
 225
 226 public:
 227   AliAnalysisTaskJetShapeHM(TString comment = "", Bool_t kMC = kFALSE);
 228  ~AliAnalysisTaskJetShapeHM();
 229
 230
 231  void SetPtJetBins(Int_t Nbin=-1, Double_t *array = 0);
 232  void SetRBins(    Int_t Nbin = 8, Double_t Rmax = 0.4) {fPsiVsRNbin = Nbin; fRmax = Rmax;}
 233  void SetPhiNbins(Int_t Nbin = 1024) {fPsiVsPhiNbin = Nbin;}
 234  void SetFilterMask(UInt_t i = 0){fFilterMask = i;}
 235  void SetEtaTrackMax(Double_t e =0.9) {fEtaTrackMax = e;}
 236  void SetPtTrackRange(Double_t min = 1., Double_t max = 100.) {fPtTrackMin = min; fPtTrackMax = max;}
 237  void SetPtJetRange(Double_t min = 1., Double_t max = 200. ) {fPtJmin  =min; fPtJmax = max; }
 238
 239  void MCprod(Bool_t kMCprod=kTRUE) {fkMCprod=kMCprod; }
 240
 241
 242  Double_t GetPtJ() {return fPtJ;}
 243  void InitHistos();
 244  Bool_t AddEvent(AliAODEvent* aodE,  AliAODJet *jet, Double_t DeltaPtJ=0.);
 245  void AddToList(TList *list);
 246
 247
 248 private:
 249  TString fComment;
 250  Bool_t fkMC;
 251  Bool_t fkMCprod;
 252  TH1F *fhEvent;//!
 253  TH1F *fhMult; //!
 254  TH1F *fhPtJ;//!
 255  TH2F *fhPtJvsPtCorr;//!
 256  TH3F *fhPsiVsR;//!
 257  TH2F *fhPsiVsRPtJ;//!
 258  TH2F *fhPhiEtaTrack;//!
 259
 260  TH3F *fhPsiVsR_MCtr; //!
 261  TH2F *fhPsiVsRPtJ_MCtr; //!
 262  TH2F *fhJetTrPtVsPartPt;//!
 263
 264
 265  TH1F *fhTMA_JAA[3];//!
 266  TH1F *fhTMA_JAp[3];//!
 267  TH1F *fhTMA_B1AA[3];//!
 268  TH1F *fhTMA_B1Ap[3];//!
 269  TH1F *fhTMA_B2AA[3];//!
 270  TH1F *fhTMA_B2Ap[3];//!
 271
 272     TH2F *fhPtresVsPt[3][2] ;//!
 273     TH2F *fhPhiresVsPhi[3][2];//!
 274     TH2F *fhEtaresVsEta[3][2];//!
 275     TH2F *fhRresVsPt[3][2];//!
 276     TH2F *fhDCAxy[3][2]; //!
 277     TH2F *fhDCAz[3][2]; //!
 278     TH3F *fhTrackPtEtaPhi[3][2];//!
 279
 280  Int_t     fPtJetNbin;
 281  TArrayD   fPtJetArray;
 282  Int_t     fPsiVsRNbin;
 283  Double_t  fRmax;
 284  Int_t     fPsiVsPhiNbin;
 285
 286  UInt_t   fFilterMask;
 287  Double_t fEtaTrackMax;
 288  Double_t fPtTrackMin;
 289  Double_t fPtTrackMax;
 290  Double_t fPtJmin;
 291  Double_t fPtJmax;
 292  Double_t fPtJ;
 293
 294  TVector3 fJvec;
 295
 296
 297  Double_t CalcR(TVector3 v1, TVector3 v2);
 298  Double_t CalcdPhi(Double_t phi1, Double_t phi2);
 299
 300  TH1F* Hist1D(const char* name, Int_t nBins, Double_t xMin, Double_t xMax,  const char* xLabel = "", Int_t color=1, const char* yLabel="");
 301  TH1F* Hist1D(const char* name, Int_t nBins, Double_t *xArray,  const char* xLabel = "", Int_t color=1, const char* yLabel="");
 302  TH2F *Hist2D(const char* name, Int_t nBinsx, Double_t xMin, Double_t xMax, Int_t nBinsy, Double_t yMin, Double_t yMax, const char* xLabel = NULL, const char* yLabel = NULL, Int_t color=1);
 303  TH2F *Hist2D(const char* name, Int_t nBinsx, Double_t xMin, Double_t xMax, Int_t nBinsy, Double_t *yArray, const char* xLabel = NULL, const char* yLabel = NULL, Int_t color=1, const char* zLabel = NULL);
 304  TH2F *Hist2D(const char* name, Int_t nBinsx, Double_t *yArrax, Int_t nBinsy, Double_t yMin, Double_t yMax, const char* xLabel = NULL, const char* yLabel = NULL, Int_t color=1, const char* zLabel = NULL);
 305   TH3F *Hist3D(const char* name, Int_t nBinsx, Double_t xMin, Double_t xMax, Int_t nBinsy, Double_t yMin, Double_t yMax, Int_t nBinsz, Double_t zMin, Double_t zMax, const char* xLabel = NULL, const char* yLabel = NULL, const char* zLabel = NULL, Int_t color=1);
 306   TH3F *Hist3D(const char* name, Int_t nBinsx, Double_t xMin, Double_t xMax, Int_t nBinsy, Double_t yMin, Double_t yMax, Int_t nBinsz, Double_t *z, const char* xLabel = NULL, const char* yLabel = NULL, const char* zLabel = NULL, Int_t color=1);
 307
 308
 309   AliAnalysisTaskJetShapeHM(const AliAnalysisTaskJetShapeHM&);            // not implemented
 310   AliAnalysisTaskJetShapeHM& operator=(const AliAnalysisTaskJetShapeHM&); // not implemented
 311
 312
 313 ClassDef(AliAnalysisTaskJetShapeHM,1)   // tbd
 314 };
 315
 316
 317
 318
 319   AliAnalysisTaskJetShape(const char *name = "");
 320   ~AliAnalysisTaskJetShape();
 321
 322   virtual void SetIsMC(Bool_t ismc=kTRUE) {fkMC=ismc;}
 323   virtual void SetCMSE(Double_t s = 7000.) {fCMSE=s;}
 324
 325
 326   virtual void     SetNonStdFile(TString c){fNonStdFile = c;}
 327   virtual void     SetBranchNames(const TString &branch1, const TString &branch2);
 328   virtual void     SetBackgroundBranch(TString &branch1, TString &branch2) { fBackgroundBranch[0] = branch1;  fBackgroundBranch[1] = branch2;};
 329
 330   virtual void     SetFilterMask(UInt_t i){fFilterMask = i;}
 331   virtual void     SetOfflineTrgMask(AliVEvent::EOfflineTriggerTypes mask) { fOfflineTrgMask = mask; }
 332   virtual void     SetCentMin(Float_t cent) { fCentMin = cent; }
 333   virtual void     SetCentMax(Float_t cent) { fCentMax = cent; }
 334   virtual void     SetEvtClassMin(Int_t evtClass) { fEvtClassMin = evtClass; }
 335   virtual void     SetEvtClassMax(Int_t evtClass) { fEvtClassMax = evtClass; }
 336   virtual void     SetJetPtCorrMin(Float_t ptJ=20, Float_t ptB=20) { fPtJcorrMin = ptJ; fPtJBcorrMin = ptB;}
 337
 338   virtual void     SetPbPb(Bool_t a = kTRUE) {fkIsPbPb = a;}
 339   // vtx
 340   virtual void     SetVtxZRange(Double_t zmin=-10., Double_t zmax=10.) {fVtxZMin = zmin; fVtxZMax = zmax;}
 341   virtual void     SetVixMinContrib(Int_t n=1) { fVtxMinContrib = n; }
 342
 343    //  virtual void     SetAngStructCloseTracks(Int_t yesno){fAngStructCloseTracks=yesno;}
 344
 345
 346   //  virtual void   ConnectInputData(Option_t *);
 347   virtual void   UserCreateOutputObjects();
 348   virtual void   UserExec(Option_t *option);
 349   virtual void   Terminate(Option_t *);
 350
 351
 352
 353
 354  private:
 355   //  TTree *fOutputTree; //!Output tree
 356
 357   TList       *fOutputList; // Output list
 358   AliESDEvent *fESD;    // ESD object
 359   AliAODEvent *fAODIn;    // AOD event
 360   AliAODEvent *fAODOut;    // AOD event
 361   AliAODExtension  *fAODExtension; //! where we take the jets from can be input or output AOD
 362
 363
 364   //  AliTriggerAnalysis * fTriggerAnalysis; // trigger analysis object, to get the offline triggers
 365
 366
 367    Bool_t   fkMC;         //is MC
 368    Double_t fCMSE;      //cms energy
 369    UInt_t fRunNb;       //run number
 370    Bool_t fkIsPhysSel;  //tbd
 371    TString       fNonStdFile; // delta AOD file
 372    UInt_t  fFilterMask;            // filter bit for slecected tracks
 373    AliVEvent::EOfflineTriggerTypes fOfflineTrgMask; // mask of offline triggers to accept
 374    Float_t fCentMin;      // lower bound on centrality
 375    Float_t fCentMax;      // upper bound on centrality
 376    Int_t   fEvtClassMin;          // lower bound on event class
 377    Int_t   fEvtClassMax;          // upper bound on event class
 378    Double_t fPtJcorrMin ;
 379    Double_t fPtJBcorrMin;
 380    Double_t fJpPtmin;
 381    Double_t fJaPtmin;
 382    Int_t fVtxMinContrib;
 383    Double_t fVtxZMin;
 384    Double_t fVtxZMax;
 385    Bool_t fkIsPbPb;
 386
 387
 388    TString fBackgroundBranch[2];  //tbd
 389    TString fJetBranchName[2]; //  name of jet branches to compare
 390
 391
 392    static const Int_t fgkbinNCent=7; // tbd
 393   Double_t fgkbinCent[fgkbinNCent+1] ;// {0, 5, 10,  20, 40, 60, 80, 100};
 394
 395   TH1F *fhPtJL ;//!
 396   TH1F *fhAreaJL ;//!
 397
 398   //  TClonesArray farray;
 399     AliAnalysisTaskJetShapeHM *fanalJetSubStrHM;//!
 400     AliAnalysisTaskJetShapeHM *fanalJetSubStrMCHM;//!
 401
 402     TH2F *fhPtresVsPt[3] ;//!
 403     TH2F *fhPhiresVsPhi[3];//!
 404     TH2F *fhEtaresVsEta[3];//!
 405     TH2F *fhDCAxy[3]; //!
 406     TH2F *fhDCAz[3]; //!
 407     TH3F *fhTrackPtEtaPhi[3];//!
 408
 409    Bool_t IsGoodEvent();
 410
 411
 412   Double_t CalcdPhi(Double_t phi1, Double_t phi2);
 413   Double_t CalcdPhi0To2pi(Double_t phi1, Double_t phi2);
 414
 415   TH1F *Hist1D(const char* name, Int_t nBins , Double_t xMin, Double_t xMax, const char* xLabel = NULL, Int_t color=1);
 416   TH2F *Hist2D(const char* name, Int_t nBinsx, Double_t xMin, Double_t xMax, Int_t nBinsy, Double_t yMin, Double_t yMax, const char* xLabel = NULL, const char* yLabel = NULL, Int_t color=1);
 417   TH3F *Hist3D(const char* name, Int_t nBinsx, Double_t xMin, Double_t xMax, Int_t nBinsy, Double_t yMin, Double_t yMax, Int_t nBinsz, Double_t zMin, Double_t zMax, const char* xLabel = NULL, const char* yLabel = NULL, const char* zLabel = NULL, Int_t color=1);
 418
 419
 420   //  AliAnalysisTaskJetShape(const AnalysisJetMain&); // not implemented
 421   //  AliAnalysisTaskJetShape& operator=(const AnalysisJetMain&); // not implemented
 422
 423
 424
 425   AliAnalysisTaskJetShape(const AliAnalysisTaskJetShape&);            // not implemented
 426   AliAnalysisTaskJetShape& operator=(const AliAnalysisTaskJetShape&); // not implemented
 427
 428
 429   ClassDef(AliAnalysisTaskJetShape, 1)
 430 };
 431
 432
 433
 434
 435
 436
 437
 438
 439
 440 #endif