Adding information from the PHOS trigger (Gustavo, Yves)
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 26 Feb 2007 16:21:03 +0000 (16:21 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 26 Feb 2007 16:21:03 +0000 (16:21 +0000)
PHOS/AliPHOSReconstructor.cxx
PHOS/AliPHOSTrigger.cxx
PHOS/AliPHOSTrigger.h
PWG4/AliTriggerPHOS.cxx
STEER/AliESD.cxx
STEER/AliESD.h
STEER/AliESDCaloCluster.cxx
STEER/AliESDCaloCluster.h

index cdd3cae..dcc2b5c 100644 (file)
@@ -34,8 +34,9 @@
 #include "AliPHOSGetter.h"
 #include "AliPHOSTracker.h"
 #include "AliRawReader.h"
+#include "AliPHOSTrigger.h"
+#include "AliPHOSGeometry.h"
 
 ClassImp(AliPHOSReconstructor)
 
 Bool_t AliPHOSReconstructor::fgDebug = kFALSE ; 
@@ -121,6 +122,62 @@ void AliPHOSReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
   
   AliDebug(2,Form("%d digits and %d rec. particles in event %d, option %s",gime->Digits()->GetEntries(),nOfRecParticles,eventNumber,GetOption()));
 
+
+  //#########Calculate trigger and set trigger info###########
+  AliPHOSTrigger tr ;
+  //   tr.SetPatchSize(1);//create 4x4 patches
+  tr.Trigger();
+  
+  Float_t maxAmp2x2  = tr.Get2x2MaxAmplitude();
+  Float_t maxAmpnxn  = tr.GetnxnMaxAmplitude();
+  Float_t ampOutOfPatch2x2  = tr.Get2x2AmpOutOfPatch() ;
+  Float_t ampOutOfPatchnxn  = tr.GetnxnAmpOutOfPatch() ;
+
+  AliPHOSGeometry * geom = gime->PHOSGeometry();
+
+  Int_t iSM2x2      = tr.Get2x2SuperModule();
+  Int_t iSMnxn      = tr.GetnxnSuperModule();
+  Int_t iCrystalPhi2x2 = tr.Get2x2CrystalPhi();
+  Int_t iCrystalPhinxn = tr.GetnxnCrystalPhi();
+  Int_t iCrystalEta2x2 = tr.Get2x2CrystalEta();
+  Int_t iCrystalEtanxn = tr.GetnxnCrystalEta();
+
+  AliDebug(2, Form("Trigger 2x2 max amp %f, out amp %f, SM %d, iphi %d ieta %d",  maxAmp2x2, ampOutOfPatch2x2, iSM2x2,iCrystalPhi2x2, iCrystalEta2x2));
+  AliDebug(2, Form("Trigger 4x4 max amp %f , out amp %f, SM %d, iphi %d, ieta %d",  maxAmpnxn, ampOutOfPatchnxn, iSMnxn,iCrystalPhinxn, iCrystalEtanxn));
+
+  Int_t iRelId2x2 []= {iSM2x2,0,iCrystalPhi2x2,iCrystalEta2x2};
+  Int_t iAbsId2x2 =-1;
+  Int_t iRelIdnxn []= {iSMnxn,0,iCrystalPhinxn,iCrystalEtanxn};
+  Int_t iAbsIdnxn =-1;
+  TVector3    pos2x2(-1,-1,-1);
+  TVector3    posnxn(-1,-1,-1);
+  geom->RelToAbsNumbering(iRelId2x2, iAbsId2x2);
+  geom->RelToAbsNumbering(iRelIdnxn, iAbsIdnxn);
+  geom->RelPosInAlice(iAbsId2x2, pos2x2);
+  geom->RelPosInAlice(iAbsIdnxn, posnxn);
+
+  TArrayF triggerPosition(6);
+  triggerPosition[0] = pos2x2(0) ;   
+  triggerPosition[1] = pos2x2(1) ;   
+  triggerPosition[2] = pos2x2(2) ;  
+  triggerPosition[3] = posnxn(0) ;   
+  triggerPosition[4] = posnxn(1) ;   
+  triggerPosition[5] = posnxn(2) ;  
+
+  TArrayF triggerAmplitudes(4);
+  triggerAmplitudes[0] = maxAmp2x2 ;   
+  triggerAmplitudes[1] = ampOutOfPatch2x2 ;    
+  triggerAmplitudes[2] = maxAmpnxn ;   
+  triggerAmplitudes[3] = ampOutOfPatchnxn ;   
+
+  //esd->SetPHOSTriggerCells(triggerPosition);
+  esd->AddPHOSTriggerPosition(triggerPosition);
+  esd->AddPHOSTriggerAmplitudes(triggerAmplitudes);
+  
+  //######################################
+
+  //Fill CaloClusters 
   for (Int_t recpart = 0 ; recpart < nOfRecParticles ; recpart++) {
     AliPHOSRecParticle * rp = dynamic_cast<AliPHOSRecParticle*>(recParticles->At(recpart));
     if (Debug()) 
@@ -130,13 +187,14 @@ void AliPHOSReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
     AliPHOSEmcRecPoint  *emcRP = gime->EmcRecPoint(ts->GetEmcIndex());
     AliESDCaloCluster   *ec    = new AliESDCaloCluster() ; 
 
-    // fills the ESDCaloCluster
+  // fills the ESDCaloCluster
     Float_t xyz[3];
     for (Int_t ixyz=0; ixyz<3; ixyz++) 
       xyz[ixyz] = rp->GetPos()[ixyz];
     
     AliDebug(2,Form("Global position xyz=(%f,%f,%f)",xyz[0],xyz[1],xyz[2]));
     
+
     Int_t  digitMult  = emcRP->GetDigitsMultiplicity();
     Int_t *digitsList = emcRP->GetDigitsList();
     UShort_t *amplList  = new UShort_t[digitMult];
@@ -156,7 +214,6 @@ void AliPHOSReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
     ec->SetClusterEnergy(rp->Energy());         //total particle energy
     ec->SetClusterDisp(emcRP->GetDispersion()); //cluster dispersion
     ec->SetPid          (rp->GetPID()) ;        //array of particle identification
-    ec->SetPrimaryIndex (rp->GetPrimaryIndex());//index of primary particle (for simulations)
     ec->SetM02(emcRP->GetM2x()) ;               //second moment M2x
     ec->SetM20(emcRP->GetM2z()) ;               //second moment M2z
     ec->SetNExMax(emcRP->GetNExMax());          //number of local maxima
@@ -168,6 +225,16 @@ void AliPHOSReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
     ec->SetClusterChi2(-1);                     //not yet implemented
     ec->SetM11(-1) ;                            //not yet implemented
 
+    //Primaries
+    ec->SetPrimaryIndex(rp->GetPrimaryIndex());
+    Int_t  primMult  = 0;
+    Int_t *primInts =  emcRP->GetPrimaries(primMult);
+    ec->SetNumberOfPrimaries(primMult);           //primary multiplicity
+    UShort_t *primList = new UShort_t[primMult];
+    for (Int_t ipr=0; ipr<primMult; ipr++) 
+      primList[ipr] = (UShort_t)(primInts[ipr]);        
+    ec->SetListOfPrimaries(primList);                  //primary List for a cluster
+    
     // add the track to the esd object
     esd->AddCaloCluster(ec);
     delete ec;    
index 2496928..f9422e3 100644 (file)
@@ -20,9 +20,9 @@
 //  Digits are grouped in TRU's (Trigger Units). A TRU consist of 16x28 
 //  crystals ordered fNTRUPhi x fNTRUZ. The algorithm searches all possible 
 //  2x2 and nxn (n multiple of 4) crystal combinations per each TRU, adding the 
-//  digits amplitude and  finding the maximum. Maxima are transformed in ADC 
-//  time samples. Each time bin is compared to the trigger threshold until it is larger 
-//  and then, triggers are set. Thresholds need to be fixed. 
+//  digits amplitude and  finding the maximum. Iti is found is maximum is isolated.
+//  Maxima are transformed in ADC time samples. Each time bin is compared to the trigger 
+//  threshold until it is larger and then, triggers are set. Thresholds need to be fixed. 
 //  Usage:
 //
 //  //Inside the event loop
 //  tr->SetL0Threshold(100);
 //  tr->SetL1JetLowPtThreshold(1000);
 //  tr->SetL1JetHighPtThreshold(20000);
+//  ....
 //  tr->Trigger(); //Execute Trigger
-//  tr->Print(""); //Print result, with "deb" option all data members 
-//  //are printed
+//  tr->Print(""); //Print data members after calculation.
+//  
 //
 //*-- Author: Gustavo Conesa & Yves Schutz (IFIC, CERN) 
 //////////////////////////////////////////////////////////////////////////////
@@ -60,7 +61,14 @@ AliPHOSTrigger::AliPHOSTrigger()
     fADCValuesHighnxn(0), fADCValuesLownxn(0),
     fADCValuesHigh2x2(0), fADCValuesLow2x2(0), fDigitsList(0),
     fL0Threshold(50), fL1JetLowPtThreshold(200), fL1JetHighPtThreshold(500),
-    fNTRU(8), fNTRUZ(2), fNTRUPhi(4),  fPatchSize(1), fSimulation(kTRUE)
+    fNTRU(8), fNTRUZ(2), fNTRUPhi(4), 
+    fNCrystalsPhi(16),
+    fNCrystalsZ(28),
+    fPatchSize(1), fIsolPatchSize(1), 
+    f2x2AmpOutOfPatch(-1), fnxnAmpOutOfPatch(-1), 
+    f2x2AmpOutOfPatchThres(2),  fnxnAmpOutOfPatchThres(2), //2 GeV out of patch 
+    fIs2x2Isol(kFALSE), fIsnxnIsol(kFALSE),  
+    fSimulation(kTRUE)
 {
   //ctor
   fADCValuesHighnxn = 0x0; //new Int_t[fTimeBins];
@@ -96,7 +104,17 @@ AliPHOSTrigger::AliPHOSTrigger(const AliPHOSTrigger & trig) :
   fNTRU(trig.fNTRU),
   fNTRUZ(trig.fNTRUZ),
   fNTRUPhi(trig.fNTRUPhi),
-  fPatchSize(trig.fPatchSize), fSimulation(trig.fSimulation)
+  fNCrystalsPhi(trig.fNCrystalsPhi),
+  fNCrystalsZ(trig. fNCrystalsZ),
+  fPatchSize(trig.fPatchSize),
+  fIsolPatchSize(trig.fIsolPatchSize), 
+  f2x2AmpOutOfPatch(trig.f2x2AmpOutOfPatch), 
+  fnxnAmpOutOfPatch(trig.fnxnAmpOutOfPatch), 
+  f2x2AmpOutOfPatchThres(trig.f2x2AmpOutOfPatchThres),  
+  fnxnAmpOutOfPatchThres(trig.fnxnAmpOutOfPatchThres), 
+  fIs2x2Isol(trig.fIs2x2Isol),
+  fIsnxnIsol(trig.fIsnxnIsol),  
+  fSimulation(trig.fSimulation)
 {
   // cpy ctor
 }
@@ -136,8 +154,6 @@ void AliPHOSTrigger::FillTRU(const TClonesArray * digits, const AliPHOSGeometry
 
   //Initilize and declare variables
   Int_t nModules     = geom->GetNModules();
-  Int_t nCrystalsPhi = geom->GetNPhi()/fNTRUPhi ;// 64/4=16
-  Int_t nCrystalsZ   = geom->GetNZ()/fNTRUZ ;// 56/2=28
   Int_t relid[4] ; 
   Float_t amp   = -1;
   Float_t timeR = -1;
@@ -145,10 +161,10 @@ void AliPHOSTrigger::FillTRU(const TClonesArray * digits, const AliPHOSGeometry
 
   //List of TRU matrices initialized to 0.
   for(Int_t k = 0; k < fNTRU*nModules ; k++){
-    TMatrixD  * amptrus   = new TMatrixD(nCrystalsPhi,nCrystalsZ) ;
-    TMatrixD  * timeRtrus = new TMatrixD(nCrystalsPhi,nCrystalsZ) ;
-    for(Int_t i = 0; i < nCrystalsPhi; i++){
-      for(Int_t j = 0; j < nCrystalsZ; j++){
+    TMatrixD  * amptrus   = new TMatrixD(fNCrystalsPhi,fNCrystalsZ) ;
+    TMatrixD  * timeRtrus = new TMatrixD(fNCrystalsPhi,fNCrystalsZ) ;
+    for(Int_t i = 0; i < fNCrystalsPhi; i++){
+      for(Int_t j = 0; j < fNCrystalsZ; j++){
        (*amptrus)(i,j)   = 0.0;
        (*timeRtrus)(i,j) = 0.0;
       }
@@ -179,12 +195,12 @@ void AliPHOSTrigger::FillTRU(const TClonesArray * digits, const AliPHOSGeometry
       //Check to which TRU in the supermodule belongs the crystal. 
       //Supermodules are divided in a TRU matrix of dimension 
       //(fNTRUPhi,fNTRUZ).
-      //Each TRU is a crystal matrix of dimension (nCrystalsPhi,nCrystalsZ)
+      //Each TRU is a crystal matrix of dimension (fNCrystalsPhi,fNCrystalsZ)
       
       //First calculate the row and column in the supermodule 
       //of the TRU to which the crystal belongs.
-      Int_t col   = (relid[3]-1)/nCrystalsZ+1; 
-      Int_t row   = (relid[2]-1)/nCrystalsPhi+1;
+      Int_t col   = (relid[3]-1)/fNCrystalsZ+1; 
+      Int_t row   = (relid[2]-1)/fNCrystalsPhi+1;
  
       //Calculate label number of the TRU  
       Int_t itru  = (row-1) + (col-1)*fNTRUPhi + (relid[0]-1)*fNTRU ;
@@ -194,8 +210,8 @@ void AliPHOSTrigger::FillTRU(const TClonesArray * digits, const AliPHOSGeometry
       TMatrixD * timeRtrus = dynamic_cast<TMatrixD *>(timeRmatrix->At(itru)) ;
 
       //Calculate row and column of the crystal inside the TRU with number itru
-      Int_t irow = (relid[2]-1) - (row-1) *  nCrystalsPhi;     
-      Int_t icol = (relid[3]-1) - (col-1) *  nCrystalsZ;
+      Int_t irow = (relid[2]-1) - (row-1) *  fNCrystalsPhi;    
+      Int_t icol = (relid[3]-1) - (col-1) *  fNCrystalsZ;
       
       (*amptrus)(irow,icol)   = amp ;
       (*timeRtrus)(irow,icol) = timeR ;
@@ -204,7 +220,7 @@ void AliPHOSTrigger::FillTRU(const TClonesArray * digits, const AliPHOSGeometry
 }
 
 //______________________________________________________________________
-void AliPHOSTrigger::GetCrystalPhiEtaIndexInModuleFromTRUIndex(const Int_t itru,const Int_t iphitru,const Int_t ietatru,Int_t &iphiMod,Int_t &ietaMod,const AliPHOSGeometry* geom) const 
+void AliPHOSTrigger::GetCrystalPhiEtaIndexInModuleFromTRUIndex(const Int_t itru,const Int_t iphitru,const Int_t ietatru,Int_t &iphiMod,Int_t &ietaMod) const 
 {
   // This method transforms the (eta,phi) index of a crystals in a 
   // TRU matrix into Super Module (eta,phi) index.
@@ -215,22 +231,89 @@ void AliPHOSTrigger::GetCrystalPhiEtaIndexInModuleFromTRUIndex(const Int_t itru,
   Int_t row = itru - (col-1)*fNTRUPhi + 1;
   
   //Calculate the (eta,phi) index in SM
-  Int_t nCrystalsPhi = geom->GetNPhi()/fNTRUPhi;
-  Int_t nCrystalsZ   = geom->GetNZ()/fNTRUZ;
   
-  iphiMod = nCrystalsPhi*(row-1) + iphitru + 1 ;
-  ietaMod = nCrystalsZ*(col-1) + ietatru + 1 ;
+  iphiMod = fNCrystalsPhi*(row-1) + iphitru + 1 ;
+  ietaMod = fNCrystalsZ*(col-1) + ietatru + 1 ;
 
 }
+
+//____________________________________________________________________________
+Bool_t AliPHOSTrigger::IsPatchIsolated(Int_t iPatchType, const TClonesArray * amptrus, const Int_t mtru, const Int_t imod, const Float_t *maxarray) {
+
+  //Calculate if the maximum patch found is isolated, find amplitude around maximum (2x2 or nxn) patch, 
+  //inside isolation patch . iPatchType = 0 means calculation for 2x2 patch, 
+  //iPatchType = 1 means calculation for nxn patch.
+  //In the next table there is an example of the different options of patch size and isolation patch size:
+  //                                                                                 Patch Size (fPatchSize)
+  //                                                             0                          1                                  2
+  //          fIsolPatchSize                 2x2 (not overlap)   4x4 (overlapped)        6x6(overlapped) ...
+  //                   1                                       4x4                      8x8                              10x10
+  //                   2                                       6x6                     12x12                           14x14    
+  //                   3                                       8x8                     16x16                           18x18
+                          
+  Bool_t b = kFALSE;
+  Float_t amp = 0;
+ //Get matrix of TRU with maximum amplitude patch.
+  Int_t itru = mtru+imod*fNTRU ; //number of tru, min 0 max 8*5.
+  TMatrixD * amptru   = dynamic_cast<TMatrixD *>(amptrus->At(itru)) ;
+  //Define patch cells
+  Int_t isolcells = fIsolPatchSize*(1+iPatchType);
+  Int_t ipatchcells = 2*(1+fPatchSize*iPatchType);
+  Int_t minrow =  static_cast<Int_t>(maxarray[1]) - isolcells;
+  Int_t mincol =  static_cast<Int_t>(maxarray[2]) - isolcells;
+  Int_t maxrow =  static_cast<Int_t>(maxarray[1]) + isolcells + ipatchcells;
+  Int_t maxcol =  static_cast<Int_t>(maxarray[2]) +  isolcells + ipatchcells;
+
+  AliDebug(2,Form("Number of added Isol Cells %d, Patch Size %d",isolcells, ipatchcells));
+  AliDebug(2,Form("Patch: minrow %d, maxrow %d, mincol %d, maxcol %d",minrow,maxrow,mincol,maxcol));
+
+  if(minrow < 0 || mincol < 0 || maxrow > fNCrystalsPhi || maxcol > fNCrystalsZ){
+    AliDebug(1,Form("Out of TRU range, cannot isolate patch"));
+    return kFALSE;
+  }
+
+  //Add amplitudes in all isolation patch
+  for(Int_t irow = minrow ; irow <  maxrow; irow ++)
+    for(Int_t icol = mincol ; icol < maxcol ; icol ++)
+      amp += (*amptru)(irow,icol);
+
+  AliDebug(2,Form("Type %d, Maximum amplitude %f, patch+isol square %f",iPatchType, maxarray[0], amp));
+
+  if(amp < maxarray[0]){
+    AliError(Form("Bad sum: Type %d, Maximum amplitude %f, patch+isol square %f",iPatchType, maxarray[0], amp));
+    return kFALSE;
+  }
+  else
+    amp-=maxarray[0]; //Calculate energy in isolation patch that do not comes from maximum patch.
+  
+  AliDebug(2, Form("Maximum amplitude %f, Out of patch %f",maxarray[0], amp));
+
+  //Fill isolation amplitude data member and say if patch is isolated.
+  if(iPatchType == 0){ //2x2 case
+    f2x2AmpOutOfPatch = amp;   
+    if(amp < f2x2AmpOutOfPatchThres)
+      b=kTRUE;
+  }
+  else  if(iPatchType == 1){ //nxn case
+    fnxnAmpOutOfPatch = amp;   
+    if(amp < fnxnAmpOutOfPatchThres)
+      b=kTRUE;
+  }
+
+  return b;
+
+}
+
+
 //____________________________________________________________________________
-void AliPHOSTrigger::MakeSlidingCell(const TClonesArray * amptrus, const TClonesArray * timeRtrus, const Int_t imod, TMatrixD *ampmax2, TMatrixD *ampmaxn, const AliPHOSGeometry *geom){
+void AliPHOSTrigger::MakeSlidingCell(const TClonesArray * amptrus, const TClonesArray * timeRtrus, const Int_t imod, TMatrixD *ampmax2, TMatrixD *ampmaxn){
   //Sums energy of all possible 2x2 (L0) and nxn (L1) crystals per each TRU. 
   //Fast signal in the experiment is given by 2x2 crystals, 
   //for this reason we loop inside the TRU crystals by 2. 
  
   //Declare and initialize varibles
-  Int_t nCrystalsPhi = geom->GetNPhi()/fNTRUPhi ;// 64/4=16
-  Int_t nCrystalsZ   = geom->GetNZ()/fNTRUZ ;// 56/2=28
   Float_t amp2 = 0 ;
   Float_t ampn = 0 ; 
   for(Int_t i = 0; i < 4; i++){
@@ -242,23 +325,22 @@ void AliPHOSTrigger::MakeSlidingCell(const TClonesArray * amptrus, const TClones
 
   //Create matrix that will contain 2x2 amplitude sums
   //used to calculate the nxn sums
-  TMatrixD  * tru2x2 = new TMatrixD(nCrystalsPhi/2,nCrystalsZ/2) ;
-  for(Int_t i = 0; i < nCrystalsPhi/2; i++)
-    for(Int_t j = 0; j < nCrystalsZ/2; j++)
-      (*tru2x2)(i,j) = 0.0;
+  TMatrixD  * tru2x2 = new TMatrixD(fNCrystalsPhi/2,fNCrystalsZ/2) ;
+  for(Int_t i = 0; i < fNCrystalsPhi/2; i++)
+    for(Int_t j = 0; j < fNCrystalsZ/2; j++)
+      (*tru2x2)(i,j) = -1.;
     
   //Loop over all TRUS in a module
-  for(Int_t itru = 0 + (imod - 1) * fNTRU ; itru < imod*fNTRU ; itru++){
+  for(Int_t itru = 0 + imod  * fNTRU ; itru < (imod+1)*fNTRU ; itru++){
     TMatrixD * amptru   = dynamic_cast<TMatrixD *>(amptrus->At(itru)) ;
     TMatrixD * timeRtru = dynamic_cast<TMatrixD *>(timeRtrus->At(itru)) ;
-    Int_t mtru = itru-(imod-1)*fNTRU ; //Number of TRU in Module
+    Int_t mtru = itru-imod*fNTRU ; //Number of TRU in Module
     
     //Sliding 2x2, add 2x2 amplitudes (NOT OVERLAP)
-    for(Int_t irow = 0 ; irow <  nCrystalsPhi; irow += 2){ 
-      for(Int_t icol = 0 ; icol < nCrystalsZ ; icol += 2){
+    for(Int_t irow = 0 ; irow <  fNCrystalsPhi; irow += 2){ 
+      for(Int_t icol = 0 ; icol < fNCrystalsZ ; icol += 2){
        amp2 = (*amptru)(irow,icol)+(*amptru)(irow+1,icol)+
          (*amptru)(irow,icol+1)+(*amptru)(irow+1,icol+1);
-       
        //Fill new matrix with added 2x2 crystals for use in nxn sums
        (*tru2x2)(irow/2,icol/2) = amp2 ;
        //Select 2x2 maximum sums to select L0 
@@ -285,10 +367,10 @@ void AliPHOSTrigger::MakeSlidingCell(const TClonesArray * amptrus, const TClones
 
     //Sliding nxn, add nxn amplitudes (OVERLAP)
     if(fPatchSize > 0){
-      for(Int_t irow = 0 ; irow <  nCrystalsPhi/2; irow++){ 
-       for(Int_t icol = 0 ; icol < nCrystalsZ/2 ; icol++){
+      for(Int_t irow = 0 ; irow <  fNCrystalsPhi/2; irow++){ 
+       for(Int_t icol = 0 ; icol < fNCrystalsZ/2 ; icol++){
          ampn = 0;
-         if( (irow+fPatchSize) < nCrystalsPhi/2 && (icol+fPatchSize) < nCrystalsZ/2){//Avoid exit the TRU
+         if( (irow+fPatchSize) < fNCrystalsPhi/2 && (icol+fPatchSize) < fNCrystalsZ/2){//Avoid exit the TRU
            for(Int_t i = 0 ; i <= fPatchSize ; i++)
              for(Int_t j = 0 ; j <= fPatchSize ; j++)
                ampn += (*tru2x2)(irow+i,icol+j);
@@ -308,7 +390,7 @@ void AliPHOSTrigger::MakeSlidingCell(const TClonesArray * amptrus, const TClones
       Int_t coln =  static_cast <Int_t> ((*ampmaxn)(2,mtru));
       for(Int_t i = 0; i<4*fPatchSize; i++){
        for(Int_t j = 0; j<4*fPatchSize; j++){
-         if( (rown+i) < nCrystalsPhi && (coln+j) < nCrystalsZ/2){//Avoid exit the TRU
+         if( (rown+i) < fNCrystalsPhi && (coln+j) < fNCrystalsZ/2){//Avoid exit the TRU
            if((*amptru)(rown+i,coln+j) > 0 &&  (*timeRtru)(rown+i,coln+j)> 0){
              if((*timeRtru)(rown+i,coln+j) <  (*ampmaxn)(3,mtru)  )
                (*ampmaxn)(3,mtru) =  (*timeRtru)(rown+i,coln+j);
@@ -325,6 +407,8 @@ void AliPHOSTrigger::MakeSlidingCell(const TClonesArray * amptrus, const TClones
       }
   }
 }
+
+
 //____________________________________________________________________________
 void AliPHOSTrigger::Print(const Option_t * opt) const 
 {
@@ -339,12 +423,15 @@ void AliPHOSTrigger::Print(const Option_t * opt) const
   printf( "               -2x2 crystals sum (not overlapped): %10.2f, in Super Module %d\n",
          f2x2MaxAmp,f2x2SM) ; 
   printf( "               -2x2 from row %d to row %d and from column %d to column %d\n", f2x2CrystalPhi, f2x2CrystalPhi+2, f2x2CrystalEta, f2x2CrystalEta+2) ; 
-
+  printf( "               -2x2 Isolation Patch %d x %d, Amplitude out of 2x2 patch is %f, threshold %f, Isolated? %d \n", 
+         2*fIsolPatchSize+2, 2*fIsolPatchSize+2,  f2x2AmpOutOfPatch,  f2x2AmpOutOfPatchThres,static_cast<Int_t> (fIs2x2Isol)) ; 
   if(fPatchSize > 0){
-    printf( "             Patch Size, n x n: %d x %d cells\n",4*fPatchSize, 4*fPatchSize);
+    printf( "             Patch Size, n x n: %d x %d cells\n",2*(fPatchSize+1), 2*(fPatchSize+1));
     printf( "               -nxn crystals sum (overlapped)    : %10.2f, in Super Module %d\n",
            fnxnMaxAmp,fnxnSM) ; 
     printf( "               -nxn from row %d to row %d and from column %d to column %d\n", fnxnCrystalPhi, fnxnCrystalPhi+4, fnxnCrystalEta, fnxnCrystalEta+4) ; 
+    printf( "               -nxn Isolation Patch %d x %d, Amplitude out of nxn patch is %f, threshold %f, Isolated? %d \n", 
+           4*fIsolPatchSize+2*(fPatchSize+1),4*fIsolPatchSize+2*(fPatchSize+1) ,  fnxnAmpOutOfPatch,  fnxnAmpOutOfPatchThres,static_cast<Int_t> (fIsnxnIsol) ) ; 
   }
   printf( "             Threshold for LO %10.1f\n", 
          fL0Threshold) ;  
@@ -367,16 +454,16 @@ void AliPHOSTrigger::Print(const Option_t * opt) const
 }
 
 //____________________________________________________________________________
-void AliPHOSTrigger::SetTriggers(const Int_t iMod, const TMatrixD * ampmax2, const TMatrixD * ampmaxn, const AliPHOSGeometry *geom)  
+void AliPHOSTrigger::SetTriggers(const TClonesArray * amptrus, const Int_t iMod, const TMatrixD * ampmax2, const TMatrixD * ampmaxn)  
 {
   //Checks the 2x2 and nxn maximum amplitude per each TRU and compares 
-  //with the different L0 and L1 triggers thresholds
+  //with the different L0 and L1 triggers thresholds. It finds if maximum amplitudes are isolated.
 
   //Initialize variables
   Float_t max2[] = {-1,-1,-1,-1} ;
   Float_t maxn[] = {-1,-1,-1,-1} ;
-  Int_t   itru2  = -1 ;
-  Int_t   itrun  = -1 ;
+  Int_t   mtru2  = -1 ;
+  Int_t   mtrun  = -1 ;
 
 
   //Find maximum summed amplitude of all the TRU 
@@ -387,14 +474,14 @@ void AliPHOSTrigger::SetTriggers(const Int_t iMod, const TMatrixD * ampmax2, con
       max2[1] =  (*ampmax2)(1,i) ; // corresponding phi position in TRU
       max2[2] =  (*ampmax2)(2,i) ; // corresponding eta position in TRU
       max2[3] =  (*ampmax2)(3,i) ; // corresponding most recent time
-      itru2   = i ; // TRU number
+      mtru2   = i ; // TRU number in module
     }
     if(maxn[0] < (*ampmaxn)(0,i) ){
       maxn[0] =  (*ampmaxn)(0,i) ; // nxn summed max amplitude
       maxn[1] =  (*ampmaxn)(1,i) ; // corresponding phi position in TRU
       maxn[2] =  (*ampmaxn)(2,i) ; // corresponding eta position in TRU
       maxn[3] =  (*ampmaxn)(3,i) ; // corresponding most recent time
-      itrun   = i ; // TRU number
+      mtrun   = i ; // TRU number in module
     }
   }
   
@@ -410,11 +497,14 @@ void AliPHOSTrigger::SetTriggers(const Int_t iMod, const TMatrixD * ampmax2, con
     f2x2MaxAmp  = max2[0] ;
     f2x2SM      = iMod ;
     maxtimeR2   = max2[3] ;
-    GetCrystalPhiEtaIndexInModuleFromTRUIndex(itru2,
+    GetCrystalPhiEtaIndexInModuleFromTRUIndex(mtru2,
                                              static_cast<Int_t>(max2[1]),
                                              static_cast<Int_t>(max2[2]),
-                                             f2x2CrystalPhi,f2x2CrystalEta,geom) ;
+                                             f2x2CrystalPhi,f2x2CrystalEta) ;
     
+    //Isolated patch?
+    fIs2x2Isol =  IsPatchIsolated(0, amptrus, mtru2, iMod, max2) ;
+
     //Transform digit amplitude in Raw Samples
     fADCValuesLow2x2  = new Int_t[nTimeBins];
     fADCValuesHigh2x2 = new Int_t[nTimeBins];
@@ -435,15 +525,18 @@ void AliPHOSTrigger::SetTriggers(const Int_t iMod, const TMatrixD * ampmax2, con
   }
 
   //Set max nxn amplitude and select L1 triggers
-  if(maxn[0] > fnxnMaxAmp ){
+  if(maxn[0] > fnxnMaxAmp  && fPatchSize > 0){
     fnxnMaxAmp  = maxn[0] ;
     fnxnSM      = iMod ;
     maxtimeRn   = maxn[3] ;
-    GetCrystalPhiEtaIndexInModuleFromTRUIndex(itrun,
+    GetCrystalPhiEtaIndexInModuleFromTRUIndex(mtrun,
                                              static_cast<Int_t>(maxn[1]),
                                              static_cast<Int_t>(maxn[2]),
-                                             fnxnCrystalPhi,fnxnCrystalEta,geom) ; 
+                                             fnxnCrystalPhi,fnxnCrystalEta) ; 
     
+    //Isolated patch?
+    fIsnxnIsol =  IsPatchIsolated(1, amptrus, mtrun, iMod, maxn) ;
+
     //Transform digit amplitude in Raw Samples
     fADCValuesHighnxn = new Int_t[nTimeBins];
     fADCValuesLownxn  = new Int_t[nTimeBins];
@@ -486,6 +579,8 @@ void AliPHOSTrigger::Trigger()
    
   //Define parameters
   Int_t nModules     = geom->GetNModules();
+  fNCrystalsPhi = geom->GetNPhi()/fNTRUPhi ;// 64/4=16
+  fNCrystalsZ   = geom->GetNZ()/fNTRUZ ;// 56/2=28
 
   //Intialize data members each time the trigger is called in event loop
   f2x2MaxAmp = -1; f2x2CrystalPhi = -1;  f2x2CrystalEta = -1;
@@ -509,10 +604,14 @@ void AliPHOSTrigger::Trigger()
   TMatrixD  * ampmax2 = new TMatrixD(4,fNTRU) ;
   TMatrixD  * ampmaxn = new TMatrixD(4,fNTRU) ;
 
-  for(Int_t imod = 1 ; imod <= nModules ; imod++) {
+  for(Int_t imod = 0 ; imod < nModules ; imod++) {
+
     //Do 2x2 and nxn sums, select maximums. 
-    MakeSlidingCell(amptrus, timeRtrus, imod, ampmax2, ampmaxn, geom);
+    MakeSlidingCell(amptrus, timeRtrus, imod, ampmax2, ampmaxn);
     //Set the trigger
-    SetTriggers(imod,ampmax2,ampmaxn, geom) ;
+    SetTriggers(amptrus,imod,ampmax2,ampmaxn) ;
   }
+
+  //Print();
+
 }
index 5b7fe03..6d5d187 100644 (file)
@@ -8,24 +8,26 @@
 
 //____________________________________________________________
 //  Class for trigger analysis.
+
+//
+//*-- Author: Gustavo Conesa & Yves Schutz (IFIC, SUBATECH, CERN)
 //  Digits are grouped in TRU's (Trigger Units). A TRU consist of 16x28 
 //  crystals ordered fNTRUPhi x fNTRUZ. The algorithm searches all possible 
-//  2x2 and nxn  (n multiple of 4) crystal combinations per each TRU, adding the 
-//  digits amplitude and finding the maximum. Maxima are transformed in ADC 
-//  time samples.  Each time bin is compared to the trigger threshold until it is larger 
-//  and then, triggers are set. Thresholds need to be fixed. 
+//  2x2 and nxn (n multiple of 4) crystal combinations per each TRU, adding the 
+//  digits amplitude and  finding the maximum. Iti is found is maximum is isolated.
+//  Maxima are transformed in ADC time samples. Each time bin is compared to the trigger 
+//  threshold until it is larger and then, triggers are set. Thresholds need to be fixed. 
 //  Usage:
 //
 //  //Inside the event loop
-//  AliPHOSTrigger *tr = new AliPHOSTrigger();//Init Trigger
+//  AliEMCALTrigger *tr = new AliEMCALTrigger();//Init Trigger
 //  tr->SetL0Threshold(100);
 //  tr->SetL1JetLowPtThreshold(1000);
 //  tr->SetL1JetHighPtThreshold(20000);
+//  ....
 //  tr->Trigger(); //Execute Trigger
-//  tr->Print("");  //Print results
-//
-//*-- Author: Gustavo Conesa & Yves Schutz (IFIC, SUBATECH, CERN)
-     
+//  tr->Print(""); //Print data members after calculation.
+//     
 // --- ROOT system ---
 
 class TClonesArray ;
@@ -64,7 +66,7 @@ class AliPHOSTrigger : public AliTriggerDetector {
   Int_t *  GetADCValuesLowGainMaxnxnSum()  {return fADCValuesLownxn; }
   Int_t *  GetADCValuesHighGainMaxnxnSum() {return fADCValuesHighnxn; }
 
-  void GetCrystalPhiEtaIndexInModuleFromTRUIndex(Int_t itru, Int_t iphitru, Int_t ietatru,Int_t &ietaMod,Int_t &iphiMod, const AliPHOSGeometry *geom) const ;
+  void GetCrystalPhiEtaIndexInModuleFromTRUIndex(Int_t itru, Int_t iphitru, Int_t ietatru,Int_t &ietaMod,Int_t &iphiMod) const ;
 
   Float_t  GetL0Threshold()          const {return fL0Threshold ; } 
   Float_t  GetL1JetLowPtThreshold()  const {return fL1JetLowPtThreshold ; }
@@ -74,7 +76,16 @@ class AliPHOSTrigger : public AliTriggerDetector {
   Int_t    GetNTRUZ()   const  {return fNTRUZ ; }
   Int_t    GetNTRUPhi() const  {return fNTRUPhi ; }
   
-  Float_t  GetPatchSize() const  {return fPatchSize ; }
+  Int_t  GetPatchSize() const  {return fPatchSize ; }
+  Int_t  GetIsolPatchSize() const  {return fIsolPatchSize ; }
+
+  Float_t  Get2x2AmpOutOfPatch() const  {return  f2x2AmpOutOfPatch; }
+  Float_t  GetnxnAmpOutOfPatch() const  {return  fnxnAmpOutOfPatch; }
+  Float_t  Get2x2AmpOutOfPatchThres() const  {return  f2x2AmpOutOfPatchThres; }
+  Float_t  GetnxnAmpOutOfPatchThres() const  {return  fnxnAmpOutOfPatchThres; }
+  Bool_t   Is2x2Isol()  const  {return  fIs2x2Isol; }
+  Bool_t   IsnxnIsol()  const  {return  fIsnxnIsol; }
+
   Bool_t   IsSimulation() const {return fSimulation ; }
 
   //Setters
@@ -95,6 +106,9 @@ class AliPHOSTrigger : public AliTriggerDetector {
     {fL1JetHighPtThreshold = amp ; }
 
   void SetPatchSize(Int_t ps)                {fPatchSize = ps ; }
+  void SetIsolPatchSize(Int_t ps)          {fIsolPatchSize = ps ; }
+  void Set2x2AmpOutOfPatchThres(Float_t th) { f2x2AmpOutOfPatchThres = th; }
+  void SetnxnAmpOutOfPatchThres(Float_t th) { fnxnAmpOutOfPatchThres = th; }
   void SetSimulation(Bool_t sim )          {fSimulation = sim ; }
 
  private:
@@ -103,9 +117,11 @@ class AliPHOSTrigger : public AliTriggerDetector {
 
   void FillTRU(const TClonesArray * digits, const AliPHOSGeometry * geom, TClonesArray * amptru, TClonesArray * timeRtru) const ;
 
-  void MakeSlidingCell(const TClonesArray * amptrus, const TClonesArray * timeRtrus, Int_t mod, TMatrixD *ampmax2, TMatrixD *ampmaxn, const AliPHOSGeometry *geom) ;
+  Bool_t IsPatchIsolated(Int_t iPatchType, const TClonesArray * amptrus, const Int_t mtru, const Int_t imod, const Float_t *maxarray) ;
+
+  void MakeSlidingCell(const TClonesArray * amptrus, const TClonesArray * timeRtrus, Int_t mod, TMatrixD *ampmax2, TMatrixD *ampmaxn) ;
 
-  void SetTriggers(Int_t iMod, const TMatrixD *ampmax2,const TMatrixD *ampmaxn, const AliPHOSGeometry *geom) ;
+ void SetTriggers(const TClonesArray * amptrus, Int_t iMod, const TMatrixD *ampmax2,const TMatrixD *ampmaxn) ;
 
  private: 
 
@@ -129,11 +145,26 @@ class AliPHOSTrigger : public AliTriggerDetector {
   Float_t fL1JetLowPtThreshold ;  //! L1 Low  pT trigger threshold
   Float_t fL1JetHighPtThreshold ; //! L1 High pT trigger threshold
 
-  Int_t   fNTRU ;                 //! Number of TRUs per module
+  Int_t   fNTRU ;                  //! Number of TRUs per module
   Int_t   fNTRUZ ;                //! Number of crystal rows per Z in one TRU
-  Int_t   fNTRUPhi ;              //! Number of crystal rows per Phi in one TRU
+  Int_t   fNTRUPhi ;             //! Number of crystal rows per Phi in one TRU
+  Int_t   fNCrystalsPhi;       //!Number of rows in a TRU
+  Int_t   fNCrystalsZ;          //!Number of columns in a TRU
+  
   Int_t fPatchSize;               //! Trigger patch factor, to be multiplied to 2x2 cells
-                                          // 0 means 2x2, 1 means nxn, 2 means 8x8 ...
+                                          // 0 means 2x2, 1 means 4x4, 2 means 6x6 ...
+  Int_t fIsolPatchSize ;      //Isolation patch size, number of rows or columns to add to 
+                                           //the 2x2 or nxn maximum amplitude patch. 
+                                           //1 means a patch around max amplitude of 2x2 of 4x4 and around         
+                                           //max ampl patch of 4x4 of 8x8 
+    
+  Float_t f2x2AmpOutOfPatch; //Amplitude in isolation cone minus maximum amplitude of the reference patch
+  Float_t fnxnAmpOutOfPatch; 
+  Float_t f2x2AmpOutOfPatchThres; //Threshold to select a trigger as isolated on f2x2AmpOutOfPatch value
+  Float_t fnxnAmpOutOfPatchThres; 
+  Float_t fIs2x2Isol; //Patch is isolated if f2x2AmpOutOfPatchThres threshold is passed
+  Float_t fIsnxnIsol ; 
+  
   Bool_t  fSimulation ;           //! Flag to do the trigger during simulation or reconstruction
   ClassDef(AliPHOSTrigger,4)
 } ;
index e69bc3c..5c08f07 100644 (file)
@@ -125,13 +125,13 @@ void AliTriggerPHOS::Exec(Option_t *)
   const Float_t aNNO   = static_cast<Float_t>(triggerAmplitudes->At(3)) ; 
 
   // trigger position
-  const TArrayI * triggerCells      = fESD->GetPHOSTriggerCells();
-  const Float_t p22    =  static_cast<Float_t>(triggerCells->At(0)) ; 
-  const Float_t phi22  =  static_cast<Float_t>(triggerCells->At(1)) ;
-  const Float_t eta22  =  static_cast<Float_t>(triggerCells->At(2)) ;
-  const Float_t pNN    =  static_cast<Float_t>(triggerCells->At(3)) ; 
-  const Float_t phiNN  =  static_cast<Float_t>(triggerCells->At(4)) ;
-  const Float_t etaNN  =  static_cast<Float_t>(triggerCells->At(5)) ; 
+  const TArrayF * triggerPosition      = fESD->GetPHOSTriggerPosition();
+  const Float_t p22    =  static_cast<Float_t>(triggerPosition->At(0)) ; 
+  const Float_t phi22  =  static_cast<Float_t>(triggerPosition->At(1)) ;
+  const Float_t eta22  =  static_cast<Float_t>(triggerPosition->At(2)) ;
+  const Float_t pNN    =  static_cast<Float_t>(triggerPosition->At(3)) ; 
+  const Float_t phiNN  =  static_cast<Float_t>(triggerPosition->At(4)) ;
+  const Float_t etaNN  =  static_cast<Float_t>(triggerPosition->At(5)) ; 
   
   Int_t       firstPhosCluster       = fESD->GetFirstPHOSCluster() ;
   const Int_t numberOfPhosClusters   = fESD->GetNumberOfPHOSClusters() ;
index 9858bab..a60b406 100644 (file)
@@ -60,10 +60,15 @@ AliESD::AliESD():
   fCaloClusters("AliESDCaloCluster",10000),
   fEMCALClusters(0), 
   fFirstEMCALCluster(-1),
+  fEMCALTriggerPosition(0x0),
+  fEMCALTriggerAmplitudes(0x0),
   fPHOSClusters(0), 
   fFirstPHOSCluster(-1),
+  fPHOSTriggerPosition(0x0),
+  fPHOSTriggerAmplitudes(0x0),
   fESDFMD(0x0),
   fESDVZERO(0x0)
+
 {
   for (Int_t i=0; i<24; i++) {
     fT0time[i] = 0;
@@ -104,8 +109,12 @@ AliESD::AliESD(const AliESD& esd):
   fCaloClusters(*((TClonesArray*)esd.fCaloClusters.Clone())),
   fEMCALClusters(esd.fEMCALClusters), 
   fFirstEMCALCluster(esd.fFirstEMCALCluster),
+  fEMCALTriggerPosition(esd. fEMCALTriggerPosition),
+  fEMCALTriggerAmplitudes(esd.fEMCALTriggerAmplitudes),
   fPHOSClusters(esd.fPHOSClusters), 
   fFirstPHOSCluster(esd.fFirstPHOSCluster),
+  fPHOSTriggerPosition(esd.fPHOSTriggerPosition),
+  fPHOSTriggerAmplitudes(esd.fPHOSTriggerAmplitudes),
   fESDFMD(esd.fESDFMD),
   fESDVZERO(esd.fESDVZERO)
 {
@@ -157,7 +166,11 @@ AliESD & AliESD::operator=(const AliESD& source) {
   fFirstPHOSCluster = source.fFirstPHOSCluster;
   fESDFMD = source.fESDFMD;
   fESDVZERO = source.fESDVZERO;
-
+  fEMCALTriggerPosition=source. fEMCALTriggerPosition;
+  fEMCALTriggerAmplitudes=source.fEMCALTriggerAmplitudes;
+  fPHOSTriggerPosition=source.fPHOSTriggerPosition;
+  fPHOSTriggerAmplitudes=source.fPHOSTriggerAmplitudes;
+  
   for (Int_t i=0; i<24; i++) {
     fT0time[i] = source.fT0time[i];
     fT0amplitude[i] = source.fT0amplitude[i];
@@ -186,6 +199,15 @@ AliESD::~AliESD()
   fCaloClusters.Delete();
   delete fESDFMD;
   delete fESDVZERO;
+//   fEMCALTriggerPosition->Delete();
+//   fEMCALTriggerAmplitudes->Delete();
+//   fPHOSTriggerPosition->Delete();
+//   fPHOSTriggerAmplitudes->Delete();
+//   delete fEMCALTriggerPosition;
+//   delete fEMCALTriggerAmplitudes;
+//   delete fPHOSTriggerPosition;
+//   delete fPHOSTriggerAmplitudes;
+
 }
 
 //______________________________________________________________________________
@@ -224,6 +246,10 @@ void AliESD::Reset()
   fPHOSClusters=0; 
   fFirstPHOSCluster=-1; 
   if (fESDFMD) fESDFMD->Clear();
+//   fEMCALTriggerPosition->Clear();
+//   fEMCALTriggerAmplitudes->Clear();
+//   fPHOSTriggerPosition->Clear();
+//   fPHOSTriggerAmplitudes->Clear();
 }
 
 Int_t AliESD::AddV0(const AliESDv0 *v) {
index cee3b1e..403da7b 100644 (file)
@@ -16,6 +16,7 @@
 
 #include <TClonesArray.h>
 #include <TObject.h>
+#include <TArrayF.h>
 
 #include "AliESDMuonTrack.h"
 #include "AliESDPmdTrack.h"
@@ -121,6 +122,11 @@ public:
     return fCaloClusters.GetEntriesFast()-1;
   }
     
+  void AddPHOSTriggerPosition(TArrayF array)   { fPHOSTriggerPosition    = new TArrayF(array) ; }
+  void AddPHOSTriggerAmplitudes(TArrayF array) { fPHOSTriggerAmplitudes  = new TArrayF(array) ; }
+  void AddEMCALTriggerPosition(TArrayF array)  { fEMCALTriggerPosition   = new TArrayF(array) ; }
+  void AddEMCALTriggerAmplitudes(TArrayF array){ fEMCALTriggerAmplitudes = new TArrayF(array) ; }
+
   void SetVertex(const AliESDVertex *vertex) {
      new (&fSPDVertex) AliESDVertex(*vertex);
   }
@@ -158,11 +164,16 @@ public:
   void  SetNumberOfEMCALClusters(Int_t clus) {fEMCALClusters = clus;}
   Int_t GetFirstEMCALCluster() const {return fFirstEMCALCluster;}
   void  SetFirstEMCALCluster(Int_t index) {fFirstEMCALCluster = index;}
+  TArrayF *GetEMCALTriggerPosition() const {return  fEMCALTriggerPosition;}
+  TArrayF *GetEMCALTriggerAmplitudes() const {return  fEMCALTriggerAmplitudes;}
 
   Int_t GetNumberOfPHOSClusters() const {return fPHOSClusters;}
   void  SetNumberOfPHOSClusters(Int_t part) { fPHOSClusters = part ; }
   void  SetFirstPHOSCluster(Int_t index) { fFirstPHOSCluster = index ; } 
   Int_t GetFirstPHOSCluster() const  { return fFirstPHOSCluster ; }
+  TArrayF *GetPHOSTriggerPosition() const {return  fPHOSTriggerPosition;}
+  TArrayF *GetPHOSTriggerAmplitudes() const {return  fPHOSTriggerAmplitudes;}
+
 
   Float_t GetT0zVertex() const {return fT0zVertex;}
   void SetT0zVertex(Float_t z) {fT0zVertex=z;}
@@ -242,10 +253,14 @@ protected:
   TClonesArray fCaloClusters;    // Calorimeter clusters for PHOS/EMCAL
   Int_t        fEMCALClusters;   // Number of EMCAL clusters (subset of caloclusters)
   Int_t        fFirstEMCALCluster; // First EMCAL cluster in the fCaloClusters list 
+  TArrayF   *fEMCALTriggerPosition; ///(x,y,z of 2x2 and x,y,z of nxn) not position of centroid but of patch corner
+  TArrayF   *fEMCALTriggerAmplitudes; //(2x2 max ampl, 2x2 amp out of patch,  nxn max ampl, nxn amp out of patch)
 
   Int_t        fPHOSClusters;     // Number of PHOS clusters (subset of caloclusters)
   Int_t        fFirstPHOSCluster; // First PHOS cluster in the fCaloClusters list 
+  TArrayF   *fPHOSTriggerPosition; //(x,y,z of 2x2 and x,y,z of nxn), not position of centroid but of patch corner
+  TArrayF   *fPHOSTriggerAmplitudes; //(2x2 max ampl, 2x2 amp out of patch,  nxn max ampl, nxn amp out of patch)
+
   AliESDFMD   *fESDFMD;   // FMD object containing rough multiplicity
   AliESDVZERO *fESDVZERO; // VZERO object containing rough multiplicity
 
index cfd2b2d..3c491c7 100644 (file)
@@ -44,6 +44,8 @@ AliESDCaloCluster::AliESDCaloCluster() :
   fM11(0),
   fNExMax(0),
   fEmcCpvDistance(9999),
+  fNumberOfPrimaries(-1),
+  fListOfPrimaries(0x0),
   fNumberOfDigits(0),
   fDigitAmplitude(0x0),
   fDigitTime(0x0),
@@ -72,6 +74,8 @@ AliESDCaloCluster::AliESDCaloCluster(const AliESDCaloCluster& clus) :
   fM11(clus.fM11),
   fNExMax(clus.fNExMax),
   fEmcCpvDistance(clus.fEmcCpvDistance),
+  fNumberOfPrimaries(clus.fNumberOfPrimaries),
+  fListOfPrimaries(0x0),
   fNumberOfDigits(clus.fNumberOfDigits),
   fDigitAmplitude(0x0),
   fDigitTime(0x0),
@@ -102,6 +106,11 @@ AliESDCaloCluster::AliESDCaloCluster(const AliESDCaloCluster& clus) :
       for (Int_t i=0; i<clus.fNumberOfDigits; i++)
        fDigitIndex[i]=clus.fDigitIndex[i];
     }
+   if (clus.fListOfPrimaries) {
+      fListOfPrimaries = new UShort_t[clus.fNumberOfPrimaries];
+      for (Int_t i=0; i<clus.fNumberOfPrimaries; i++)
+       fListOfPrimaries[i]=clus.fListOfPrimaries[i];
+    }
   }
 }
 
@@ -125,6 +134,10 @@ AliESDCaloCluster &AliESDCaloCluster::operator=(const AliESDCaloCluster& source)
   fM11 = source.fM11;
   fNExMax = source.fNExMax;
   fEmcCpvDistance = source.fEmcCpvDistance;
+
+  fNumberOfPrimaries = source.fNumberOfPrimaries;
+  delete fListOfPrimaries; fListOfPrimaries=0x0;
+
   fNumberOfDigits = source.fNumberOfDigits;
   delete fDigitAmplitude; fDigitAmplitude=0x0;
   delete fDigitTime; fDigitTime = 0x0;
@@ -152,6 +165,11 @@ AliESDCaloCluster &AliESDCaloCluster::operator=(const AliESDCaloCluster& source)
       for (Int_t i=0; i<source.fNumberOfDigits; i++)
        fDigitIndex[i]=source.fDigitIndex[i];
     }
+   if (source.fListOfPrimaries) {
+      fListOfPrimaries = new UShort_t[source.fNumberOfPrimaries];
+      for (Int_t i=0; i<source.fNumberOfPrimaries; i++)
+       fListOfPrimaries[i]=source.fListOfPrimaries[i];
+    }
   }
 
   return *this;
@@ -166,7 +184,7 @@ AliESDCaloCluster::~AliESDCaloCluster(){
   //
   // AliESDCaloCluster is the owner of the arrays
   // even if they are created outside
-
+  delete[] fListOfPrimaries;
   delete[] fDigitAmplitude;
   delete[] fDigitTime;
   delete[] fDigitIndex;
index 1ba3bed..787253c 100644 (file)
@@ -81,6 +81,12 @@ public:
   void SetEmcCpvDistance(Float_t dEmcCpv) { fEmcCpvDistance = dEmcCpv; }
   Float_t GetEmcCpvDistance() const       { return fEmcCpvDistance; }
 
+  void SetNumberOfPrimaries(Int_t nprim)      { fNumberOfPrimaries = nprim; }
+  Int_t GetNumberOfPrimaries() const         { return fNumberOfPrimaries; }
+  
+  void SetListOfPrimaries(UShort_t *prim)       { fListOfPrimaries = prim;}
+  UShort_t *GetListOfPrimaries() const          { return fListOfPrimaries;}
+  
   void SetNumberOfDigits(Int_t ndig)      { fNumberOfDigits = ndig; }
   Int_t GetNumberOfDigits() const         { return fNumberOfDigits; }
   
@@ -113,7 +119,8 @@ protected:
   UShort_t  fNExMax ;          // number of (Ex-)maxima before unfolding
   Float_t   fEmcCpvDistance;   // the distance from PHOS EMC rec.point to the closest CPV rec.point
 
-
+  Int_t     fNumberOfPrimaries;   // number of primaries that generated the cluster
+  UShort_t*   fListOfPrimaries;   //[ fNumberOfPrimaries] list of primaries that generated the cluster
 
   Int_t     fNumberOfDigits;   // number of calorimeter digits in cluster
                                // Very important! The streamer needs to