MUON/AliMUONDigitizer.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2000, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   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  **************************************************************************/
  15
  16 /*
  17 $Log$
  18 Revision 1.8  2002/02/13 09:03:24  jchudoba
  19 Remove some deletes from dtor, those objects are deleted earlier in Exec() method (where they are created)
  20
  21 Revision 1.7  2001/11/22 11:15:41  jchudoba
  22 Proper deletion of arrays (thanks to Rene Brun)
  23
  24 Revision 1.6  2001/11/02 12:55:45  jchudoba
  25 cleanup of the code, add const to Get methods
  26
  27 Revision 1.4  2001/10/18 14:44:09  jchudoba
  28 Define constant MAXTRACKS for maximum number of tracks associated with 1 digit
  29
  30 Revision 1.3  2001/10/04 20:01:54  jchudoba
  31 changes for TTask implementation, some other small editing
  32
  33 Revision 1.2  2001/07/28 10:46:04  hristov
  34 AliRunDigitizer.h included; typos corrected
  35
  36 Revision 1.1  2001/07/27 15:41:01  jchudoba
  37 merging/digitization classes
  38
  39 */
  40 #include <TTree.h>
  41 #include <TObjArray.h>
  42 #include <TFile.h>
  43 #include <TDirectory.h>
  44 #include <iostream.h>
  45
  46 #include "AliMUONDigitizer.h"
  47 #include "AliMUONConstants.h"
  48 #include "AliMUONChamber.h"
  49 #include "AliMUONHitMapA1.h"
  50 #include "AliMUON.h"
  51 #include "AliMUONHit.h"
  52 #include "AliMUONPadHit.h"
  53 #include "AliMUONDigit.h"
  54 #include "AliMUONTransientDigit.h"
  55 #include "AliRun.h"
  56 #include "AliPDG.h"
  57 #include "AliRunDigitizer.h"
  58
  59 ClassImp(AliMUONDigitizer)
  60
  61 //___________________________________________
  62 AliMUONDigitizer::AliMUONDigitizer() :AliDigitizer()
  63 {
  64 // Default ctor - don't use it
  65   fHits = 0;
  66   fPadHits = 0;
  67   fHitMap = 0;
  68   fTDList = 0;
  69 }
  70
  71 //___________________________________________
  72 AliMUONDigitizer::AliMUONDigitizer(AliRunDigitizer* manager)
  73     :AliDigitizer(manager)
  74 {
  75 // ctor which should be used
  76   fHitMap  = 0;
  77   fTDList  = 0;
  78   fHits    = 0;
  79   fPadHits = 0;
  80   fDebug   = 0;
  81   if (GetDebug()>2)
  82     cerr<<"AliMUONDigitizer::AliMUONDigitizer"
  83         <<"(AliRunDigitizer* manager) was processed"<<endl;
  84 }
  85
  86 //------------------------------------------------------------------------
  87 AliMUONDigitizer::~AliMUONDigitizer()
  88 {
  89 // Destructor
  90   if (fHits)       delete fHits;
  91   if (fPadHits)    delete fPadHits;
  92 }
  93
  94 //------------------------------------------------------------------------
  95 Bool_t AliMUONDigitizer::Exists(const AliMUONPadHit *padhit) const
  96 {
  97   return (fHitMap[fNch]->TestHit(padhit->PadX(),padhit->PadY()));
  98 }
  99
 100 //------------------------------------------------------------------------
 101 void AliMUONDigitizer::Update(AliMUONPadHit *padhit)
 102 {
 103     AliMUONTransientDigit *pdigit =
 104       static_cast<AliMUONTransientDigit*>(
 105       fHitMap[fNch]->GetHit(padhit->PadX(),padhit->PadY()));
 106
 107     // update charge
 108     //
 109     Int_t iqpad    = padhit->QPad();        // charge per pad
 110     pdigit->AddSignal(iqpad);
 111     pdigit->AddPhysicsSignal(iqpad);
 112
 113     // update list of tracks
 114     //
 115     Int_t track, charge;
 116     track = fTrack+fMask;
 117     if (fSignal) {
 118       charge = iqpad;
 119     } else {
 120       charge = kBgTag;
 121     }
 122     pdigit->UpdateTrackList(track,charge);
 123 }
 124
 125 //------------------------------------------------------------------------
 126 void AliMUONDigitizer::CreateNew(AliMUONPadHit *padhit)
 127 {
 128 // Create new AliMUONTransientDigit and add it to the fTDList
 129
 130   fTDList->AddAtAndExpand(
 131     new AliMUONTransientDigit(fNch,fDigits),fCounter);
 132   fHitMap[fNch]->SetHit(padhit->PadX(),padhit->PadY(),fCounter);
 133   AliMUONTransientDigit* pdigit =
 134     (AliMUONTransientDigit*)fTDList->Last();
 135   // list of tracks
 136   Int_t track, charge;
 137   track = fTrack+fMask;
 138   if (fSignal) {
 139     charge = padhit->QPad();
 140   } else {
 141     charge = kBgTag;
 142   }
 143   pdigit->AddToTrackList(track,charge);
 144   fCounter++;
 145 }
 146
 147
 148 //------------------------------------------------------------------------
 149 Bool_t AliMUONDigitizer::Init()
 150 {
 151 // Initialization
 152
 153   fHits     = new TClonesArray("AliMUONHit",1000);
 154   fPadHits  = new TClonesArray("AliMUONPadHit",1000);
 155   return kTRUE;
 156 }
 157
 158
 159 //------------------------------------------------------------------------
 160 //void AliMUONDigitizer::Digitize()
 161 void AliMUONDigitizer::Exec(Option_t* option)
 162 {
 163
 164   TString optionString = option;
 165   if (optionString.Data() == "deb") {
 166     cout<<"AliMUONDigitizer::Exec: called with option deb "<<endl;
 167     fDebug = 3;
 168   }
 169   AliMUONChamber*   iChamber;
 170   AliSegmentation*  segmentation;
 171
 172   if (GetDebug()>2) cerr<<"   AliMUONDigitizer::Digitize() starts"<<endl;
 173   fTDList = new TObjArray;
 174
 175   AliMUON *pMUON  = (AliMUON *) gAlice->GetModule("MUON");
 176   if (!pMUON) {
 177     cerr<<"AliMUONDigitizer::Digitize Error:"
 178         <<" module MUON not found in the input file"<<endl;
 179     return;
 180   }
 181   pMUON->MakeBranchInTreeD(fManager->GetTreeD());
 182   fHitMap= new AliMUONHitMapA1* [AliMUONConstants::NCh()];
 183
 184   //
 185   // loop over cathodes
 186   //
 187
 188   for (int icat = 0; icat < 2; icat++) {
 189     fCounter = 0;
 190     for (Int_t i = 0; i < AliMUONConstants::NCh(); i++) {
 191       iChamber = &(pMUON->Chamber(i));
 192 //      if (!(iChamber->Nsec() == 1 && icat == 1)) {
 193         segmentation = iChamber->SegmentationModel(icat+1);
 194         fHitMap[i] = new AliMUONHitMapA1(segmentation, fTDList);
 195 //      }
 196     }
 197
 198
 199 // Loop over files to digitize
 200     fSignal = kTRUE;
 201     for (Int_t inputFile=0; inputFile<fManager->GetNinputs();
 202          inputFile++) {
 203 // Connect MUON branches
 204
 205       if (inputFile > 0 ) fSignal = kFALSE;
 206       TBranch *branchHits = 0;
 207       TBranch *branchPadHits = 0;
 208       TTree *treeH = fManager->GetInputTreeH(inputFile);
 209       if (GetDebug()>2) {
 210         cerr<<"   inputFile , cathode = "<<inputFile<<" "
 211             <<icat<<endl;
 212         cerr<<"   treeH, fHits "<<treeH<<" "<<fHits<<endl;
 213       }
 214       if (treeH && fHits) {
 215         branchHits = treeH->GetBranch("MUON");
 216         if (branchHits) {
 217           fHits->Clear();
 218           branchHits->SetAddress(&fHits);
 219         }
 220         else
 221           Error("Exec","branch MUON was not found");
 222       }
 223       if (GetDebug()>2) cerr<<"   branchHits = "<<branchHits<<endl;
 224
 225       if (treeH && fPadHits) {
 226         branchPadHits = treeH->GetBranch("MUONCluster");
 227         if (branchPadHits)
 228           branchPadHits->SetAddress(&fPadHits);
 229         else
 230           Error("Exec","branch MUONCluster was not found");
 231       }
 232       if (GetDebug()>2) cerr<<"   branchPadHits = "<<branchPadHits<<endl;
 233
 234 //
 235 //   Loop over tracks
 236 //
 237
 238       Int_t ntracks = (Int_t) treeH->GetEntries();
 239
 240       for (fTrack = 0; fTrack < ntracks; fTrack++) {
 241         if (GetDebug()>2) cerr<<"   fTrack = "<<fTrack<<endl;
 242         fHits->Clear();
 243         fPadHits->Clear();
 244         branchHits->GetEntry(fTrack);
 245         branchPadHits->GetEntry(fTrack);
 246
 247 //
 248 //   Loop over hits
 249
 250         AliMUONHit* mHit;
 251         for(Int_t i = 0; i < fHits->GetEntriesFast(); ++i) {
 252           mHit = static_cast<AliMUONHit*>(fHits->At(i));
 253           fNch = mHit->Chamber()-1;  // chamber number
 254           if (fNch > AliMUONConstants::NCh()-1) {
 255             cerr<<"AliMUONDigitizer: ERROR: "
 256                 <<"fNch > AliMUONConstants::NCh()-1, fNch, NCh(): "
 257                 <<fNch<<", "<< AliMUONConstants::NCh()<<endl;
 258             return;
 259           }
 260           iChamber = &(pMUON->Chamber(fNch));
 261 //
 262 // Loop over pad hits
 263           for (AliMUONPadHit* mPad =
 264                  (AliMUONPadHit*)pMUON->FirstPad(mHit,fPadHits);
 265                mPad;
 266                mPad = (AliMUONPadHit*)pMUON->NextPad(fPadHits))
 267           {
 268             Int_t cathode  = mPad->Cathode();      // cathode number
 269             Int_t ipx      = mPad->PadX();         // pad number on X
 270             Int_t ipy      = mPad->PadY();         // pad number on Y
 271             Int_t iqpad    = Int_t(mPad->QPad());  // charge per pad
 272             if (cathode != (icat+1)) continue;
 273
 274             fMask = fManager->GetMask(inputFile);
 275             fDigits[0] = ipx;
 276             fDigits[1] = ipy;
 277             if (!(fHitMap[fNch]->ValidateHit(fDigits[0], fDigits[1]))) continue;
 278             fDigits[2] = icat;
 279             fDigits[3] = iqpad;
 280             if (inputFile == 0) {
 281               fDigits[4] = iqpad;
 282             } else {
 283               fDigits[4] = 0;
 284             }
 285             if (mHit->Particle() == kMuonPlus ||
 286                 mHit->Particle() == kMuonMinus) {
 287               fDigits[5] = (mPad->HitNumber()) + fMask;
 288             } else fDigits[5] = -1;
 289
 290             // build the list of fired pads and update the info,
 291             // fDigits is input for Update(mPad)
 292
 293             if (!Exists(mPad)) {
 294               CreateNew(mPad);
 295             } else {
 296               Update(mPad);
 297             } //  end if Exists(mPad)
 298           } //end loop over clusters
 299         } // hit loop
 300       } // track loop
 301     } // end file loop
 302     if (GetDebug()>2) cerr<<"END OF FILE LOOP"<<endl;
 303
 304     Int_t tracks[kMAXTRACKS];
 305     Int_t charges[kMAXTRACKS];
 306     Int_t nentries = fTDList->GetEntriesFast();
 307
 308     for (Int_t nent = 0; nent < nentries; nent++) {
 309       AliMUONTransientDigit *address = (AliMUONTransientDigit*)fTDList->At(nent);
 310       if (address == 0) continue;
 311       Int_t ich = address->Chamber();
 312       Int_t   q = address->Signal();
 313       iChamber = &(pMUON->Chamber(ich));
 314 //
 315 //  Digit Response (noise, threshold, saturation, ...)
 316       AliMUONResponse * response = iChamber->ResponseModel();
 317       q = response->DigitResponse(q);
 318
 319       if (!q) continue;
 320
 321       fDigits[0] = address->PadX();
 322       fDigits[1] = address->PadY();
 323       fDigits[2] = address->Cathode();
 324       fDigits[3] = q;
 325       fDigits[4] = address->Physics();
 326       fDigits[5] = address->Hit();
 327
 328       Int_t nptracks = address->GetNTracks();
 329
 330       if (nptracks > kMAXTRACKS) {
 331         if (GetDebug() >0) {
 332           cerr<<"AliMUONDigitizer: nptracks > 10 "<<nptracks;
 333           cerr<<"reset to max value "<<kMAXTRACKS<<endl;
 334         }
 335         nptracks = kMAXTRACKS;
 336       }
 337       if (nptracks > 2 && GetDebug() >2) {
 338         cerr<<"AliMUONDigitizer: nptracks > 2 "<<nptracks<<endl;
 339         printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,fDigits[0],fDigits[1],q);
 340       }
 341       for (Int_t tr = 0; tr < nptracks; tr++) {
 342         tracks[tr]   = address->GetTrack(tr);
 343         charges[tr]  = address->GetCharge(tr);
 344       }      //end loop over list of tracks for one pad
 345       // Sort list of tracks according to charge
 346       if (nptracks > 1) {
 347         SortTracks(tracks,charges,nptracks);
 348       }
 349       if (nptracks < kMAXTRACKS ) {
 350         for (Int_t i = nptracks; i < kMAXTRACKS; i++) {
 351           tracks[i]  = 0;
 352           charges[i] = 0;
 353         }
 354       }
 355
 356       // fill digits
 357       pMUON->AddDigits(ich,tracks,charges,fDigits);
 358     }
 359
 360     fManager->GetTreeD()->Fill();
 361
 362     pMUON->ResetDigits();  //
 363     fTDList->Clear();
 364
 365
 366     for(Int_t ii = 0; ii < AliMUONConstants::NCh(); ++ii) {
 367       if (fHitMap[ii]) {
 368         delete fHitMap[ii];
 369         fHitMap[ii] = 0;
 370       }
 371     }
 372   } //end loop over cathodes
 373   if (GetDebug()>2)
 374     cerr<<"AliMUONDigitizer::Exec: writing the TreeD: "
 375         <<fManager->GetTreeD()->GetName()<<endl;
 376   fManager->GetTreeD()->Write(0,TObject::kOverwrite);
 377   delete [] fHitMap;
 378   delete fTDList;
 379
 380   if (fHits)    fHits->Delete();
 381   if (fPadHits) fPadHits->Delete();
 382 }
 383
 384
 385 //------------------------------------------------------------------------
 386 void AliMUONDigitizer::SortTracks(Int_t *tracks,Int_t *charges,Int_t ntr)
 387 {
 388   //
 389   // Sort the list of tracks contributing to a given digit
 390   // Only the 3 most significant tracks are acctually sorted
 391   //
 392
 393   //
 394   //  Loop over signals, only 3 times
 395   //
 396
 397   Int_t qmax;
 398   Int_t jmax;
 399   Int_t idx[3] = {-2,-2,-2};
 400   Int_t jch[3] = {-2,-2,-2};
 401   Int_t jtr[3] = {-2,-2,-2};
 402   Int_t i,j,imax;
 403
 404   if (ntr<3) imax=ntr;
 405   else imax=3;
 406   for(i=0;i<imax;i++){
 407     qmax=0;
 408     jmax=0;
 409
 410     for(j=0;j<ntr;j++){
 411
 412       if((i == 1 && j == idx[i-1])
 413          ||(i == 2 && (j == idx[i-1] || j == idx[i-2]))) continue;
 414
 415       if(charges[j] > qmax) {
 416         qmax = charges[j];
 417         jmax=j;
 418       }
 419     }
 420
 421     if(qmax > 0) {
 422       idx[i]=jmax;
 423       jch[i]=charges[jmax];
 424       jtr[i]=tracks[jmax];
 425     }
 426
 427   }
 428
 429   for(i=0;i<3;i++){
 430     if (jtr[i] == -2) {
 431       charges[i]=0;
 432       tracks[i]=0;
 433     } else {
 434       charges[i]=jch[i];
 435       tracks[i]=jtr[i];
 436     }
 437   }
 438 }