PWG/EMCAL/AliEmcalTriggerMaker.cxx

   1 // $Id$
   2 //
   3 // Class to make emcal particles in AOD/ESD events.
   4 //
   5 // Author: J.Kral
   6
   7 #include <TClonesArray.h>
   8 #include <iostream>
   9
  10 #include "AliLog.h"
  11 #include "AliEmcalTriggerPatchInfo.h"
  12 #include "AliEmcalTriggerSetupInfo.h"
  13
  14 #include "AliEmcalTriggerMaker.h"
  15 #include "AliEMCALTriggerTypes.h"
  16 #include "AliEMCALGeometry.h"
  17
  18 #include "AliVCaloTrigger.h"
  19 #include "AliVCaloCells.h"
  20
  21 ClassImp(AliEmcalTriggerMaker)
  22
  23 using namespace std;
  24
  25 //________________________________________________________________________
  26 AliEmcalTriggerMaker::AliEmcalTriggerMaker() :
  27   AliAnalysisTaskEmcal("AliEmcalTriggerMaker",kFALSE),
  28   fCaloTriggersOutName("EmcalTriggers"),
  29   fCaloTriggerSetupOutName("EmcalTriggersSetup"),
  30   fCaloTriggersOut(0),
  31   fCaloTriggerSetupOut(0)
  32 {
  33   // Constructor.
  34 }
  35
  36 //________________________________________________________________________
  37 AliEmcalTriggerMaker::AliEmcalTriggerMaker(const char *name) :
  38   AliAnalysisTaskEmcal(name,kFALSE),
  39   fCaloTriggersOutName("EmcalTriggers"),
  40   fCaloTriggerSetupOutName("EmcalTriggersSetup"),
  41   fCaloTriggersOut(0),
  42   fCaloTriggerSetupOut(0)
  43 {
  44   // Constructor.
  45
  46 }
  47
  48 //________________________________________________________________________
  49 AliEmcalTriggerMaker::~AliEmcalTriggerMaker()
  50 {
  51   // Destructor.
  52 }
  53
  54 //________________________________________________________________________
  55 void AliEmcalTriggerMaker::ExecOnce()
  56 {
  57   // Init the analysis.
  58
  59   AliAnalysisTaskEmcal::ExecOnce();
  60
  61   if (!fInitialized)
  62     return;
  63
  64   if (!fCaloTriggersOutName.IsNull()) {
  65     fCaloTriggersOut = new TClonesArray("AliEmcalTriggerPatchInfo");
  66     fCaloTriggersOut->SetName(fCaloTriggersOutName);
  67
  68     if (!(InputEvent()->FindListObject(fCaloTriggersOutName))) {
  69       InputEvent()->AddObject(fCaloTriggersOut);
  70     }
  71     else {
  72       fInitialized = kFALSE;
  73       AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fCaloTriggersOutName.Data()));
  74       return;
  75     }
  76   }
  77
  78   if (!fCaloTriggerSetupOutName.IsNull()) {
  79     fCaloTriggerSetupOut = new AliEmcalTriggerSetupInfo();
  80     fCaloTriggerSetupOut->SetName(fCaloTriggerSetupOutName);
  81
  82     if (!(InputEvent()->FindListObject(fCaloTriggerSetupOutName))) {
  83       InputEvent()->AddObject(fCaloTriggerSetupOut);
  84     }
  85     else {
  86       fInitialized = kFALSE;
  87       AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fCaloTriggerSetupOutName.Data()));
  88       return;
  89     }
  90   }
  91
  92 }
  93
  94 //________________________________________________________________________
  95 Bool_t AliEmcalTriggerMaker::Run()
  96 {
  97   // Create the emcal particles
  98
  99   Int_t globCol, globRow, tBits, cellAbsId[4];
 100   Int_t absId, adcAmp;
 101   Int_t i, j, k, iMain, cmCol, cmRow, cmiCellCol, cmiCellRow;
 102   Int_t jetTrigger, iTriggers;
 103   Double_t amp, ca, eMain, cmiCol, cmiRow;
 104
 105   TVector3 centerGeo, center1, center2, centerMass, edge1, edge2, vertex;
 106
 107   AliEmcalTriggerPatchInfo *trigger;
 108
 109   // delete patch array, clear setup object
 110   fCaloTriggersOut->Delete();
 111   fCaloTriggerSetupOut->Clean();
 112
 113   if( !fCaloTriggers ){
 114     AliError(Form("Calo triggers container %s not available.", fCaloTriggersName.Data()));
 115     return kTRUE;
 116   }
 117   if( !fCaloCells ){
 118     AliError(Form("Calo cells container %s not available.", fCaloCellsName.Data()));
 119     return kTRUE;
 120   }
 121
 122   // must reset before usage, or the class will fail
 123   fCaloTriggers->Reset();
 124
 125   // dig out common data (thresholds)
 126   // 0 - jet high, 1 - gamma high, 2 - jet low, 3 - gamma low
 127   fCaloTriggerSetupOut->SetThresholds( fCaloTriggers->GetL1Threshold( 0 ),
 128                                        fCaloTriggers->GetL1Threshold( 1 ),
 129                                        fCaloTriggers->GetL1Threshold( 2 ),
 130                                        fCaloTriggers->GetL1Threshold( 3 ));
 131
 132   // class is not empty
 133   if( fCaloTriggers->GetEntries() > 0 ){
 134
 135     iTriggers = 0;
 136     iMain = -1;
 137     eMain = -1;
 138
 139     // save primary vertex in vector
 140     vertex.SetXYZ( fVertex[0], fVertex[1], fVertex[2] );
 141
 142     // go throuth the trigger channels
 143     while( fCaloTriggers->Next() ){
 144
 145       // check if jet trigger low or high
 146       fCaloTriggers->GetTriggerBits( tBits );
 147
 148       jetTrigger = 0;
 149       if(( tBits >> ( kTriggerTypeEnd + kL1JetLow )) & 1 )
 150         jetTrigger = 1;
 151       if(( tBits >> ( kTriggerTypeEnd + kL1JetHigh )) & 1)
 152         jetTrigger = jetTrigger | 2;
 153
 154       if( jetTrigger == 0 )
 155         continue;
 156
 157       // get position in global 2x2 tower coordinates
 158       // A0 left bottom (0,0)
 159       fCaloTriggers->GetPosition( globCol, globRow );
 160
 161       // get the absolute trigger ID
 162       fGeom->GetAbsFastORIndexFromPositionInEMCAL( globCol, globRow, absId );
 163       // convert to the 4 absId of the cells composing the trigger channel
 164       fGeom->GetCellIndexFromFastORIndex( absId, cellAbsId );
 165
 166       // get low left edge (eta max, phi min)
 167       fGeom->GetGlobal( cellAbsId[0], edge1 );
 168
 169       // sum the available energy in the 32/32 window of cells
 170       // step over trigger channels and get all the corresponding cells
 171       // make CM
 172       amp = 0;
 173       cmiCol = 0;
 174       cmiRow = 0;
 175       for( i = 0; i < 16; i++ ){
 176         for( j = 0; j < 16; j++ ){
 177           // get the 4 cells composing the trigger channel
 178           fGeom->GetAbsFastORIndexFromPositionInEMCAL( globCol+i, globRow+j, absId );
 179           fGeom->GetCellIndexFromFastORIndex( absId, cellAbsId );
 180           // add amplitudes and find patch edges
 181           for( k = 0; k < 4; k++ ){
 182             ca = fCaloCells->GetCellAmplitude( cellAbsId[k] );
 183             amp += ca;
 184             cmiCol += ca*(Double_t)i;
 185             cmiRow += ca*(Double_t)j;
 186           }
 187         }
 188       } // 32x32 cell window
 189       if( amp == 0 ){
 190         AliDebug(2,"EMCal trigger patch with 0 energy.");
 191         continue;
 192       }
 193
 194       // get the CM and patch index
 195       cmiCol /= amp;
 196       cmiRow /= amp;
 197       cmCol = globCol + (Int_t)cmiCol;
 198       cmRow = globRow + (Int_t)cmiRow;
 199
 200       // get the patch and corresponding cells
 201       fGeom->GetAbsFastORIndexFromPositionInEMCAL( cmCol, cmRow, absId );
 202       fGeom->GetCellIndexFromFastORIndex( absId, cellAbsId );
 203
 204       // find which out of the 4 cells is closest to CM and get it's position
 205       cmiCellCol = TMath::Nint( cmiCol * 2. );
 206       cmiCellRow = TMath::Nint( cmiRow * 2. );
 207       fGeom->GetGlobal( cellAbsId[(cmiCellRow%2)*2 + cmiCellCol%2], centerMass );
 208
 209       // get up right edge (eta min, phi max)
 210       // get the absolute trigger ID
 211       fGeom->GetAbsFastORIndexFromPositionInEMCAL( globCol+15, globRow+15, absId );
 212       // convert to the 4 absId of the cells composing the trigger channel
 213       fGeom->GetCellIndexFromFastORIndex( absId, cellAbsId );
 214
 215       fGeom->GetGlobal( cellAbsId[3], edge2 );
 216
 217       // get the geometrical center as an average of two diagonally
 218       // adjacent patches in the center
 219       // picking two diagonally closest cells from the patches
 220       fGeom->GetAbsFastORIndexFromPositionInEMCAL( globCol+7, globRow+7, absId );
 221       fGeom->GetCellIndexFromFastORIndex( absId, cellAbsId );
 222       fGeom->GetGlobal( cellAbsId[3], center1 );
 223
 224       fGeom->GetAbsFastORIndexFromPositionInEMCAL( globCol+8, globRow+8, absId );
 225       fGeom->GetCellIndexFromFastORIndex( absId, cellAbsId );
 226       fGeom->GetGlobal( cellAbsId[0], center2 );
 227
 228       centerGeo = center1;
 229       centerGeo += center2;
 230       centerGeo *= 0.5;
 231
 232       // relate all to primary vertex
 233       centerGeo -= vertex;
 234       centerMass -= vertex;
 235       edge1 -= vertex;
 236       edge2 -= vertex;
 237
 238       // get ADC amplitude
 239       fCaloTriggers->GetL1TimeSum( adcAmp );
 240
 241       // save the trigger object
 242       new ((*fCaloTriggersOut)[iTriggers])AliEmcalTriggerPatchInfo();
 243       trigger = (AliEmcalTriggerPatchInfo*)fCaloTriggersOut->At( iTriggers );
 244       iTriggers++;
 245
 246       trigger->SetCenterGeo( centerGeo, amp );
 247       trigger->SetCenterMass( centerMass, amp );
 248       trigger->SetEdge1( edge1, amp );
 249       trigger->SetEdge2( edge2, amp );
 250       trigger->SetADCAmp( adcAmp );
 251       trigger->SetTriggerBits( tBits );
 252
 253       // check if more energetic than others for main patch marking
 254       if( eMain < amp ){
 255         eMain = amp;
 256         iMain = iTriggers - 1;
 257       }
 258
 259 //       cout << " pi:" << trigger->GetPhiMin() << " px:" << trigger->GetPhiMax();
 260 //       cout << " pg:" << trigger->GetPhiGeo() << " " << (trigger->GetPhiMin()+trigger->GetPhiMax()) / 2.;
 261 //       cout << " pc:" << trigger->GetPhiCM();
 262 //       cout << " ei:" << trigger->GetEtaMin() << " ex:" << trigger->GetEtaMax();
 263 //       cout << " eg:" << trigger->GetEtaGeo() << " " << (trigger->GetEtaMin()+trigger->GetEtaMax()) / 2.;
 264 //       cout << " ec:" << trigger->GetEtaCM();
 265 //       cout << " e:" << trigger->GetPatchE();
 266 //       cout << " jl:" << trigger->IsJetLow() << " jh:" << trigger->IsJetHigh() << endl;
 267
 268     } // triggers
 269
 270     // mark the most energetic patch as main
 271     if( iMain > -1 ){
 272       trigger = (AliEmcalTriggerPatchInfo*)fCaloTriggersOut->At( iMain );
 273       tBits = trigger->GetTriggerBits();
 274       // main trigger flag
 275       tBits = tBits | ( 1 << 24 );
 276       trigger->SetTriggerBits( tBits );
 277     }
 278
 279   } // there are some triggers
 280
 281   return kTRUE;
 282 }