JETAN/jetan2004/AliTkConeJet.cxx

   1 //$Id$
   2
   3 #include <Riostream.h>
   4 #include <TROOT.h>
   5 #include <TClonesArray.h>
   6 #include <TParticle.h>
   7 #include <TMath.h>
   8
   9 #include "AliTkTowerV2.h"
  10 #include "AliTkConeJet.h"
  11
  12 #define Ntowers__ 100
  13
  14 AliTkConeJet::AliTkConeJet() : TObject(),
  15                          fEta(-999),fPhi(-999),fEt(-999),fType(0),
  16                          fNTowers(0),fTowers(0),fPtCut(0.),
  17                          fCEta(-999),fCPhi(-999),fCEt(-999),
  18                          fPLength(0.),fXAxis(0.),fYAxis(0.),fZAxis(0.),
  19                          fPtLead(0.),fLeadPart(0),fNParts(0),fParts(0)
  20 {
  21   fTowers = new TClonesArray("AliTkTowerV2",Ntowers__);
  22 }
  23
  24 AliTkConeJet::AliTkConeJet(Float_t pt, Float_t eta, Float_t phi,Int_t type)
  25                         : TObject(),
  26                           fEta(eta),fPhi(phi),fEt(pt),fType(type),
  27                           fNTowers(0),fTowers(0),fPtCut(0.),
  28                           fCEta(eta),fCPhi(phi),fCEt(pt),
  29                           fPLength(0.),fXAxis(0.),fYAxis(0.),fZAxis(0.),
  30                           fPtLead(0.),fLeadPart(0),fNParts(0),fParts(0)
  31 {
  32   fTowers = new TClonesArray("AliTkTowerV2",0);
  33 }
  34
  35 AliTkConeJet::AliTkConeJet(const AliTkConeJet &j)
  36                         : TObject(),
  37                           fEta(-999), fPhi(-999), fEt(-999),fType(0),
  38                           fNTowers(0),fTowers(0),fPtCut(0.),
  39                           fCEta(-999),fCPhi(-999),fCEt(-999),
  40                           fPLength(0.),fXAxis(0.),fYAxis(0.),fZAxis(0.),
  41                           fPtLead(0.),fLeadPart(0),fNParts(0),fParts(0)
  42 {
  43   fEta = j.getEta();
  44   fPhi = j.getPhi();
  45   fEt = j.getEt();
  46   fCEta = j.getCEta();
  47   fCPhi = j.getCPhi();
  48   fCEt = j.getCEt();
  49   fNTowers = j.getNTowers();
  50   fPtCut = j.getPtCut();
  51   fCEta = j.getCEta();
  52   fCPhi = j.getCPhi();
  53   fCEt = j.getCEt();
  54   fPLength = j.getPLength();
  55   fXAxis = j.getXAxis();
  56   fYAxis = j.getYAxis();
  57   fZAxis = j.getZAxis();
  58   fPtLead = j.getPtLead();
  59   fLeadPart = new TParticle(*j.getLeadPart());
  60   fNParts = j.getNParts();
  61
  62   // need to build a copy of the old TClonesArray...
  63   fTowers = new TClonesArray("AliTkTowerV2",Ntowers__);
  64   TClonesArray *otherTowers = j.getTowerList();
  65   AliTkTowerV2 *myTower;
  66   TIterator *iter = otherTowers->MakeIterator();
  67   Int_t i =0;
  68   while ((myTower = (AliTkTowerV2 *) iter->Next()) != NULL) {
  69     new ((*fTowers)[i]) AliTkTowerV2(*myTower);
  70     i++;
  71   }
  72   if(i!=fNTowers)
  73     cerr << "AliTkConeJet: should not happen " << i << " " << fNTowers << endl;
  74   delete iter;
  75
  76   // need to build a copy of the old TClonesArray...
  77   fParts = new TClonesArray("TParticle",j.getNParts());
  78   TClonesArray *otherParts = j.getParts();
  79   TParticle *myParticle;
  80   iter = otherParts->MakeIterator();
  81   i =0;
  82   while ((myParticle = (TParticle *) iter->Next()) != NULL) {
  83     new ((*fParts)[i]) TParticle(*myParticle);
  84     i++;
  85   }
  86   if(i!=fNParts)
  87     cerr << "AliTkConeJet: should not happen " << i << " " << fNParts << endl;
  88   delete iter;
  89 }
  90
  91 AliTkConeJet::~AliTkConeJet()
  92 {
  93   delete fTowers;
  94   if(fLeadPart) delete fLeadPart;
  95   if(fParts) delete fParts;
  96 }
  97
  98 void AliTkConeJet::addTower(AliTkTowerV2 *tower)
  99 {
 100   new ((*fTowers)[fNTowers]) AliTkTowerV2(*tower);
 101   fNTowers++;
 102 }
 103
 104 void AliTkConeJet::Clear(Option_t *option)
 105 {
 106   TObject::Clear(option);
 107   fTowers->Delete();
 108   fNTowers=0;
 109   if(fParts) delete fParts;
 110   fParts=0;
 111   if(fLeadPart) delete fLeadPart;
 112   fLeadPart=0;
 113   fEta=-999;
 114   fPhi=-999;
 115   fEt=-999;
 116   fPtCut=0.;
 117   fCEta=-999;
 118   fCPhi=-999;
 119   fCEt=-999;
 120   fPLength=0.;
 121   fXAxis=0.;
 122   fYAxis=0.;
 123   fZAxis=0.;
 124   fPtLead=0.;
 125   fLeadPart=0;
 126   fNParts=0;
 127   fParts=0;
 128 }
 129
 130 void AliTkConeJet::calculateFromParticles(Float_t &et, Float_t &eta, Float_t &phi, Float_t ptcut)
 131 {
 132   Float_t px=0.,py=0.,pz=0.;
 133   TIterator *titer = fTowers->MakeIterator();
 134   AliTkTowerV2 *tower = NULL;
 135   while ((tower = (AliTkTowerV2 *)titer->Next()) != NULL) {
 136     TClonesArray *particles = tower->getParticleList();
 137     TParticle *particle = NULL;
 138     TIterator *part = particles->MakeIterator();
 139     while ((particle = (TParticle *)part->Next()) != NULL) {
 140       if(particle->Pt()<ptcut) continue;
 141       px+=particle->Px();
 142       py+=particle->Py();
 143       pz+=particle->Pz();
 144     }
 145     delete part;
 146   }
 147   delete titer;
 148
 149   et = TMath::Sqrt(px*px+py*py);
 150   Float_t p = TMath::Sqrt(px*px+py*py+pz*pz);
 151   Float_t theta = (pz==0)?TMath::PiOver2():TMath::ACos(pz/p);
 152   phi = TMath::Pi()+TMath::ATan2(-py,-px);
 153   eta = -TMath::Log(TMath::Tan(theta/2.));
 154 }
 155
 156 /* or take e = p_t * cosh(eta) */
 157 Float_t AliTkConeJet::getE(Float_t ptcut) const
 158 {
 159   Float_t e=0.;
 160
 161   TIterator *iter = fTowers->MakeIterator();
 162   AliTkTowerV2 *tower = NULL;
 163   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 164     TClonesArray *particles = tower->getParticleList();
 165     TParticle *particle = NULL;
 166     TIterator *part = particles->MakeIterator();
 167     while ((particle = (TParticle *)part->Next()) != NULL) {
 168       if(particle->Pt()<ptcut) continue;
 169       e+=particle->P();
 170     }
 171     delete part;
 172   }
 173   delete iter;
 174   return e;
 175 }
 176
 177 /* or take e = p_t * cosh(eta) */
 178 Float_t AliTkConeJet::getECharged(Float_t ptcut) const
 179 {
 180   Float_t e=0.;
 181
 182   TIterator *iter = fTowers->MakeIterator();
 183   AliTkTowerV2 *tower = NULL;
 184   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 185     TClonesArray *particles = tower->getChargedParticleList();
 186     TParticle *particle = NULL;
 187     TIterator *part = particles->MakeIterator();
 188     while ((particle = (TParticle *)part->Next()) != NULL) {
 189       if(particle->Pt()<ptcut) continue;
 190       e+=particle->P();
 191     }
 192     particles->Delete();
 193     delete particles;
 194     delete part;
 195   }
 196   delete iter;
 197   return e;
 198 }
 199
 200 Float_t AliTkConeJet::getEtCharged() const
 201 {
 202   Float_t EtCharged = 0;
 203   TIterator *iter = fTowers->MakeIterator();
 204   AliTkTowerV2 *tower = NULL;
 205   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 206     EtCharged += tower->getEtCharged();
 207   }
 208   delete iter;
 209   return EtCharged;
 210 }
 211
 212 /* use to get fake rate if
 213    in mixed case particles were marked */
 214 Float_t AliTkConeJet::getEtChargedMarked(Float_t ptcut) const
 215 {
 216   Float_t et=0.;
 217
 218   TIterator *iter = fTowers->MakeIterator();
 219   AliTkTowerV2 *tower = NULL;
 220   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 221     TClonesArray *particles = tower->getChargedParticleList();
 222     TParticle *particle = NULL;
 223     TIterator *part = particles->MakeIterator();
 224     while ((particle = (TParticle *)part->Next()) != NULL) {
 225       if(particle->GetWeight()>-100) continue; //(-123)
 226       if(particle->Pt()<ptcut) continue;
 227       et+=particle->Pt();
 228     }
 229     particles->Delete();
 230     delete particles;
 231    delete part;
 232   }
 233   delete iter;
 234   return et;
 235 }
 236
 237 /* use to get fake rate if
 238    in mixed case particles were marked */
 239 Float_t AliTkConeJet::getEtMarked(Float_t ptcut) const
 240 {
 241   Float_t et=0.;
 242   TIterator *iter = fTowers->MakeIterator();
 243   AliTkTowerV2 *tower = NULL;
 244   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 245     TClonesArray *particles = tower->getParticleList();
 246     TParticle *particle = NULL;
 247     TIterator *part = particles->MakeIterator();
 248     while ((particle = (TParticle *)part->Next()) != NULL) {
 249       if(particle->GetWeight()>-100) continue; //(-123)
 250       if(particle->Pt()<ptcut) continue;
 251       et+=particle->Pt();
 252     }
 253     delete part;
 254   }
 255   delete iter;
 256   return et;
 257 }
 258
 259 Float_t AliTkConeJet::getEtMarkedFrac(Float_t ptcut) const
 260 {
 261   Float_t et=0.,etall=0.;
 262   TIterator *iter = fTowers->MakeIterator();
 263   AliTkTowerV2 *tower = NULL;
 264   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 265     TClonesArray *particles = tower->getParticleList();
 266     TParticle *particle = NULL;
 267     TIterator *part = particles->MakeIterator();
 268     while ((particle = (TParticle *)part->Next()) != NULL) {
 269       if(particle->Pt()<ptcut) continue;
 270       etall+=particle->Pt();
 271       if(particle->GetWeight()>-100) continue; //(-123)
 272       et+=particle->Pt();
 273     }
 274     delete part;
 275   }
 276   delete iter;
 277   if(etall>0) return et/etall;
 278   return 1;
 279 }
 280
 281 Float_t AliTkConeJet::getEtEM() const
 282 {
 283   Float_t EtEM = 0;
 284   TIterator *iter = fTowers->MakeIterator();
 285   AliTkTowerV2 *tower = NULL;
 286   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 287     EtEM += tower->getEtEM();
 288   }
 289   delete iter;
 290   return EtEM;
 291 }
 292
 293 TClonesArray *AliTkConeJet::getParticles(Float_t ptcut) const
 294 {
 295   TClonesArray *allParticles = new TClonesArray("TParticle",10000);
 296   Int_t nParticles = 0;
 297   TIterator *iter = fTowers->MakeIterator();
 298   AliTkTowerV2 *tower = NULL;
 299   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 300     TClonesArray *particles = tower->getParticleList();
 301     TParticle *particle = NULL;
 302     TIterator *part = particles->MakeIterator();
 303     while ((particle = (TParticle *)part->Next()) != NULL) {
 304       if(particle->Pt()<ptcut) continue;
 305       new ((*allParticles)[nParticles]) TParticle(*particle);
 306       nParticles++;
 307     }
 308     delete part;
 309   }
 310   delete iter;
 311   return allParticles;
 312 }
 313
 314 TClonesArray *AliTkConeJet::getChargedParticles(Float_t ptcut) const
 315 {
 316   TClonesArray *chargedParticles = new TClonesArray("TParticle",10000);
 317   Int_t nChargedParticles = 0;
 318   TIterator *iter = fTowers->MakeIterator();
 319   AliTkTowerV2 *tower = NULL;
 320   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 321     TClonesArray *particles = tower->getChargedParticleList();
 322     TParticle *particle = NULL;
 323     TIterator *part = particles->MakeIterator();
 324     while ((particle = (TParticle *)part->Next()) != NULL) {
 325       if(particle->Pt()<ptcut) continue;
 326       new ((*chargedParticles)[nChargedParticles]) TParticle(*particle);
 327       nChargedParticles++;
 328     }
 329     particles->Delete();
 330     delete particles;
 331     delete part;
 332   }
 333   delete iter;
 334   return chargedParticles;
 335 }
 336
 337 TClonesArray *AliTkConeJet::getNeutralParticles(Float_t ptcut) const
 338 {
 339   TClonesArray *neutralParticles = new TClonesArray("TParticle",10000);
 340   Int_t nNeutralParticles = 0;
 341   TIterator *iter = fTowers->MakeIterator();
 342   AliTkTowerV2 *tower = NULL;
 343   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 344     TClonesArray *particles = tower->getNeutralParticleList();
 345     TParticle *particle = NULL;
 346     TIterator *part = particles->MakeIterator();
 347     while ((particle = (TParticle *)part->Next()) != NULL) {
 348       if(particle->Pt()<ptcut) continue;
 349       new ((*neutralParticles)[nNeutralParticles]) TParticle(*particle);
 350       nNeutralParticles++;
 351     }
 352     particles->Delete();
 353     delete particles;
 354     delete part;
 355   }
 356   delete iter;
 357   return neutralParticles;
 358 }
 359
 360 Int_t AliTkConeJet::getNParticles() const
 361 {
 362   Int_t n = -1;
 363   TClonesArray *p = getParticles();
 364   n = p->GetEntries();
 365   p->Delete();
 366   delete p;
 367   return n;
 368 }
 369
 370 Int_t AliTkConeJet::getNChargedParticles() const
 371 {
 372   Int_t n = -1;
 373   TClonesArray *p = getChargedParticles();
 374   n = p->GetEntries();
 375   p->Delete();
 376   delete p;
 377   return n;
 378 }
 379
 380 Int_t AliTkConeJet::getNNeutralParticles() const
 381 {
 382   Int_t n = -1;
 383   TClonesArray *p = getNeutralParticles();
 384   n = p->GetEntries();
 385   p->Delete();
 386   delete p;
 387   return n;
 388 }
 389
 390 Int_t AliTkConeJet::getNParticles(Float_t ptCut) const
 391 {
 392   Int_t n = 0;
 393   TClonesArray *particles = getParticles();
 394   TIterator *iter = particles->MakeIterator();
 395   TParticle *particle;
 396   while ((particle = (TParticle *)iter->Next()) != NULL) {
 397     if (particle->Pt() > ptCut) {
 398       n++;
 399     }
 400   }
 401   delete iter;
 402   particles->Delete();
 403   delete particles;
 404   return n;
 405 }
 406
 407 Int_t AliTkConeJet::getNChargedParticles(Float_t ptCut) const
 408 {
 409   Int_t n = 0;
 410   TClonesArray *particles = getChargedParticles();
 411   TIterator *iter = particles->MakeIterator();
 412   TParticle *particle;
 413   while ((particle = (TParticle *)iter->Next()) != NULL) {
 414     if (particle->Pt() > ptCut) {
 415       n++;
 416     }
 417   }
 418   delete iter;
 419   particles->Delete();
 420   delete particles;
 421   return n;
 422 }
 423
 424 Int_t AliTkConeJet::getNNeutralParticles(Float_t ptCut) const
 425 {
 426   Int_t n = 0;
 427   TClonesArray *particles = getNeutralParticles();
 428   TIterator *iter = particles->MakeIterator();
 429   TParticle *particle;
 430   while ((particle = (TParticle *)iter->Next()) != NULL) {
 431     if (particle->Pt() > ptCut) {
 432       n++;
 433     }
 434   }
 435   delete iter;
 436   particles->Delete();
 437   delete particles;
 438   return n;
 439 }
 440
 441 void AliTkConeJet::getAxis(Float_t &x,Float_t &y,Float_t &z,Float_t ptcut) const
 442 {
 443   x=0.;
 444   y=0.;
 445   z=0.;
 446
 447   TIterator *iter = fTowers->MakeIterator();
 448   AliTkTowerV2 *tower = NULL;
 449   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 450     TClonesArray *particles = tower->getParticleList();
 451     TParticle *particle = NULL;
 452     TIterator *part = particles->MakeIterator();
 453     while ((particle = (TParticle *)part->Next()) != NULL) {
 454       if(particle->Pt()<ptcut) continue;
 455       x+=particle->Px();
 456       y+=particle->Py();
 457       z+=particle->Pz();
 458     }
 459     delete part;
 460   }
 461   delete iter;
 462 }
 463
 464 void AliTkConeJet::getChAxis(Float_t &x,Float_t &y,Float_t &z,Float_t ptcut) const
 465 {
 466   x=0.;
 467   y=0.;
 468   z=0.;
 469
 470   TIterator *iter = fTowers->MakeIterator();
 471   AliTkTowerV2 *tower = NULL;
 472   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 473     TClonesArray *particles = tower->getChargedParticleList();
 474     TParticle *particle = NULL;
 475     TIterator *part = particles->MakeIterator();
 476     while ((particle = (TParticle *)part->Next()) != NULL) {
 477       if(particle->Pt()<ptcut) continue;
 478       x+=particle->Px();
 479       y+=particle->Py();
 480       z+=particle->Pz();
 481     }
 482     particles->Delete();
 483     delete particles;
 484     delete part;
 485   }
 486   delete iter;
 487   Float_t length=TMath::Sqrt(x*x+y*y+z*z);
 488   if(length!=0){
 489     x/=length;
 490     y/=length;
 491     z/=length;
 492   }
 493 }
 494
 495 TParticle* AliTkConeJet::getLeadingPart(Float_t ptcut) const
 496 {
 497   TParticle* plead=new TParticle;
 498   Float_t ptlead=0;
 499
 500   TIterator *iter = fTowers->MakeIterator();
 501   AliTkTowerV2 *tower = NULL;
 502   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 503     TClonesArray *particles = tower->getParticleList();
 504     TParticle *particle = NULL;
 505     TIterator *part = particles->MakeIterator();
 506     while ((particle = (TParticle *)part->Next()) != NULL) {
 507       if(particle->Pt()<ptcut) continue;
 508       if(particle->Pt()>ptlead){
 509         ptlead=particle->Pt();
 510         delete plead;
 511         plead=new TParticle(*particle);
 512       }
 513     }
 514     delete part;
 515   }
 516   delete iter;
 517   return plead;
 518 }
 519
 520 TParticle* AliTkConeJet::getLeadingChPart(Float_t ptcut) const
 521 {
 522   TParticle* plead=new TParticle;
 523   Float_t ptlead=0;
 524
 525   TIterator *iter = fTowers->MakeIterator();
 526   AliTkTowerV2 *tower = NULL;
 527   while ((tower = (AliTkTowerV2 *)iter->Next()) != NULL) {
 528     TClonesArray *particles = tower->getChargedParticleList();
 529     TParticle *particle = NULL;
 530     TIterator *part = particles->MakeIterator();
 531     while ((particle = (TParticle *)part->Next()) != NULL) {
 532       if(particle->Pt()<ptcut) continue;
 533       if(particle->Pt()>ptlead){
 534         ptlead=particle->Pt();
 535         delete plead;
 536         plead=new TParticle(*particle);
 537       }
 538     }
 539     particles->Delete();
 540     delete particles;
 541     delete part;
 542   }
 543   delete iter;
 544   return plead;
 545 }
 546
 547 void AliTkConeJet::calculateValues(Float_t ptcut)
 548 {
 549   fPtCut=ptcut;
 550   if(fParts) delete fParts;
 551   fParts=new TClonesArray("TParticle",fNTowers*25);
 552
 553   Float_t px=0.,py=0.,pz=0.;
 554   TParticle* plead=new TParticle;
 555   Float_t ptlead=0;
 556   Int_t counter=0;
 557   TIterator *titer = fTowers->MakeIterator();
 558   AliTkTowerV2 *tower = NULL;
 559   while ((tower = (AliTkTowerV2 *)titer->Next()) != NULL) {
 560     TClonesArray *particles = tower->getParticleList();
 561
 562     TParticle *particle = NULL;
 563     TIterator *part = particles->MakeIterator();
 564     while ((particle = (TParticle *)part->Next()) != NULL) {
 565       if(particle->Pt()<ptcut) continue;
 566       px+=particle->Px();
 567       py+=particle->Py();
 568       pz+=particle->Pz();
 569       new ((*fParts)[counter]) TParticle(*particle);
 570       counter++;
 571       if(particle->Pt()>ptlead){
 572         ptlead=particle->Pt();
 573         delete plead;
 574         plead=new TParticle(*particle);
 575       }
 576     }
 577     delete part;
 578   }
 579   delete titer;
 580
 581   fCEt = TMath::Sqrt(px*px+py*py);
 582   Float_t p = TMath::Sqrt(px*px+py*py+pz*pz);
 583   fPLength = p;
 584   //Float_t theta = ((pz==0)||(p==0))?TMath::PiOver2():TMath::ACos(pz/p);
 585   fCPhi = TMath::Pi()+TMath::ATan2(-py,-px);
 586   //fCEta = (theta==0)?-100:TMath::Log(TMath::Tan(theta/2.));
 587   fCEta = (p==pz)?-100:0.5*TMath::Log((p+pz)/(p-pz));
 588   if(p){
 589     fXAxis= px/p;
 590     fYAxis= py/p;
 591     fZAxis= pz/p;
 592   }
 593   fPtLead=ptlead;
 594   if(fLeadPart) delete fLeadPart;
 595   fLeadPart=new TParticle(*plead);
 596   delete plead;
 597   fNParts=counter;
 598 }
 599
 600 Int_t AliTkConeJet::Compare(const TObject *obj) const
 601 {
 602   Double_t val=((AliTkConeJet*)obj)->getEt();
 603
 604   if(fEt>val) return -1; //qsort is ascending
 605   else if (fEt<val) return 1;
 606   else return 0;
 607 }
 608
 609 ostream& operator<<(ostream& s,const AliTkConeJet &p)
 610 {
 611   return s << "AliTkConeJet info: eta=" << p.getEta()
 612            << " phi=" << p.getPhi() << " Et=" << p.getEt();
 613 }
 614
 615 Float_t AliTkConeJet::Diff(const AliTkConeJet *jet1, const AliTkConeJet *jet2, Float_t &etdiff, Float_t &phidiff, Float_t &etadiff)
 616 {
 617   Float_t ret=0;
 618   phidiff=TMath::Abs(jet1->getPhi()-jet2->getPhi());
 619   if(phidiff>TMath::Pi()) ret=TMath::TwoPi()-phidiff;
 620   etadiff=TMath::Abs(jet1->getEta()-jet2->getEta());
 621   ret=TMath::Sqrt(phidiff*phidiff+etadiff*etadiff);
 622   etdiff=jet1->getEt()-jet2->getEt();
 623   return ret;
 624 }
 625
 626 ClassImp(AliTkConeJet)