Updated SDigitizer; Added AliTOFanalyzeSDigits.C macro
authorvicinanz <vicinanz@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 3 May 2002 07:34:19 +0000 (07:34 +0000)
committervicinanz <vicinanz@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 3 May 2002 07:34:19 +0000 (07:34 +0000)
TOF/AliTOF.cxx
TOF/AliTOF.h
TOF/AliTOFSDigitizer.cxx
TOF/AliTOFSDigitizer.h
TOF/AliTOFanalyzeSDigits.C
TOF/AliTOFhits2sdigits.C

index bb9d23e..9d2c324 100644 (file)
@@ -15,6 +15,9 @@
 
 /*
 $Log$
+Revision 1.36  2002/04/19 14:40:51  vicinanz
+Updated SDigitizer
+
 Revision 1.35  2002/03/21 13:52:53  vicinanz
 Minor changes to AliTOF constructor
 
@@ -1002,3 +1005,20 @@ void AliTOF::Raw2Digits(Int_t evNumber)
   tD->Write(0,TObject::kOverwrite);
 } 
 
+////////////////////////////////////////////////////////////////////////
+void AliTOF::RecreateSDigitsArray() {
+//
+// delete TClonesArray fSDigits and create it again
+//  needed for backward compatability with PPR test production
+//
+  delete fSDigits;
+  fSDigits       = new TClonesArray("AliTOFSDigit",  1000);
+}
+////////////////////////////////////////////////////////////////////////
+void AliTOF::CreateSDigitsArray() {
+//
+// create TClonesArray fSDigits
+//  needed for backward compatability with PPR test production
+//
+  fSDigits       = new TClonesArray("AliTOFSDigit",  1000);
+}
index 7ec09bc..36311c4 100644 (file)
@@ -83,13 +83,16 @@ public:
   virtual void    ResetDigits();
   TClonesArray *SDigits() const {return fSDigits;}
   TClonesArray *ReconParticles() const {return fReconParticles;}
+  void RecreateSDigitsArray();
+  void CreateSDigitsArray();
+
   Int_t   fNevents ;        // Number of events to digitize
 
 protected:
   TFolder* fFGeom ;       //  Folder that holds the Geometry definition
   TTask*   fDTask ;       //  TOF Digitizer container
   TTask*   fReTask;       //  TOF Reconstructioner container
-  TClonesArray* fSDigits; // List of summable digits
+  TClonesArray* fSDigits; //! List of summable digits
   Int_t    fNSDigits;           // Number of sdigits
   TClonesArray* fReconParticles; // List of reconstructed particles
   AliTOFMerger *fMerger;   // ! pointer to merger
@@ -126,7 +129,7 @@ protected:
 
 private:
 
-  ClassDef(AliTOF,3)  // Time Of Flight base class
+  ClassDef(AliTOF,4)  // Time Of Flight base class
 };
  
 #endif /* ALITOF_H */
index 7bc9fdc..66d1198 100644 (file)
@@ -61,6 +61,8 @@ ClassImp(AliTOFSDigitizer)
   fEvent1=0;
   fEvent2=0;
   ftail    = 0;
+  fSelectedSector=0;
+  fSelectedPlate =0;
 }
            
 //____________________________________________________________________________ 
@@ -69,6 +71,8 @@ ClassImp(AliTOFSDigitizer)
   fEvent1=evNumber1;
   fEvent2=fEvent1+nEvents;
   ftail    = 0;
+  fSelectedSector=0; // by default we sdigitize all sectors
+  fSelectedPlate =0; // by default we sdigitize all plates in all sectors
 
   fHeadersFile = HeaderFile ; // input filename (with hits)
   TFile * file = (TFile*) gROOT->GetFile(fHeadersFile.Data() ) ;
@@ -164,6 +168,17 @@ void AliTOFSDigitizer::Exec(Option_t *verboseOption, Option_t *allEvents) {
     return;
   }
 
+  // is pointer to fSDigits non zero after changes?
+  cout<<"TOF fSDigits pointer:"<<TOF->SDigits()<<endl;
+
+  // recreate TClonesArray fSDigits - for backward compatibility
+  if (TOF->SDigits() == 0) {
+    TOF->CreateSDigitsArray();
+  } else {
+    TOF->RecreateSDigitsArray();
+  }
+
+
   if (fEdgeTails) ftail = new TF1("tail",TimeWithTail,-2,2,3);
 
   Int_t nselectedHits=0;
@@ -180,6 +195,8 @@ void AliTOFSDigitizer::Exec(Option_t *verboseOption, Option_t *allEvents) {
     fEvent2= (Int_t) gAlice->TreeE()->GetEntries();
   }
 
+  Bool_t thereIsNotASelection=(fSelectedSector==0) && (fSelectedPlate==0);
+
   for (Int_t ievent = fEvent1; ievent < fEvent2; ievent++) {
     cout << "------------------- "<< GetName() << " -------------" << endl ;
     cout << "Sdigitizing event " << ievent << endl;
@@ -215,7 +232,7 @@ void AliTOFSDigitizer::Exec(Option_t *verboseOption, Option_t *allEvents) {
     // create hit map
     AliTOFHitMap *hitMap = new AliTOFHitMap(TOF->SDigits());
 
-    // decrease required CPU time
+    // increase performances in terms of CPU time
     TH->SetBranchStatus("*",0); // switch off all branches
     TH->SetBranchStatus("TOF*",1); // switch on only TOF
 
@@ -242,96 +259,101 @@ void AliTOFSDigitizer::Exec(Option_t *verboseOption, Option_t *allEvents) {
        Int_t tracknum = tofHit->GetTrack();
        vol[0] = tofHit->GetSector();
        vol[1] = tofHit->GetPlate();
-       vol[2] = tofHit->GetStrip();
-       vol[3] = tofHit->GetPadx();
-       vol[4] = tofHit->GetPadz();
-
-       Bool_t dummy=((tracknum==previousTrack) && (vol[0]==previousSector) && (vol[1]==previousPlate) && (vol[2]==previousStrip));
-
-       Bool_t isCloneOfThePrevious=dummy && ((vol[3]==previousPadX) && (vol[4]==previousPadZ));
-
-       // much stronger check to be inserted here
-
-       if(!isCloneOfThePrevious){
-         // update "previous" values
-         // in fact, we are yet in the future, so the present is past
-         previousTrack=tracknum;
-         previousSector=vol[0];
-         previousPlate=vol[1];
-         previousStrip=vol[2];
-         previousPadX=vol[3];
-         previousPadZ=vol[4];
-
-         nselectedHits++;
-         nselectedHitsinEv++;
-         if (particle->GetFirstMother() < 0){
-           nHitsFromPrim++;
-         } // counts hits due to primary particles
-
-         Float_t Xpad = tofHit->GetDx();
-         Float_t Zpad = tofHit->GetDz();
-         Float_t xStrip=AliTOFConstants::fgkXPad*(vol[3]-0.5-0.5*AliTOFConstants::fgkNpadX)+Xpad;
-         Float_t zStrip=AliTOFConstants::fgkZPad*(vol[4]-0.5-0.5*AliTOFConstants::fgkNpadZ)+Zpad;
-         Float_t geantTime = tofHit->GetTof(); // unit [s]
-         geantTime *= 1.e+09;  // conversion from [s] to [ns]
-
-         //cout << "geantTime " << geantTime << " [ns]" << endl;
-         Int_t nActivatedPads = 0, nFiredPads = 0;
-         Bool_t isFired[4] = {kFALSE, kFALSE, kFALSE, kFALSE};
-         Float_t tofAfterSimul[4] = {0., 0., 0., 0.};
-         Float_t qInduced[4] = {0.,0.,0.,0.};
-         Int_t nPlace[4] = {0, 0, 0, 0};
-         Float_t averageTime = 0.;
-         SimulateDetectorResponse(zStrip,xStrip,geantTime,nActivatedPads,nFiredPads,isFired,nPlace,qInduced,tofAfterSimul,averageTime);
-         if(nFiredPads) {
-           for(Int_t indexOfPad=0; indexOfPad<nActivatedPads; indexOfPad++) {
-             if(isFired[indexOfPad]){ // the pad has fired
-               Float_t timediff=geantTime-tofAfterSimul[indexOfPad];
-
-               if(timediff>=0.2) nlargeTofDiff++;
-
-               digit[0] = (Int_t) ((tofAfterSimul[indexOfPad]*1.e+03)/fTdcBin); // TDC bin number (each bin -> 25. ps)
-
-               Float_t landauFactor = gRandom->Landau(fAdcMean, fAdcRms); 
-               digit[1] = (Int_t) (qInduced[indexOfPad] * landauFactor); // ADC bins (each bin -> 0.25 (or 0.03) pC)
-
-               // recalculate the volume only for neighbouring pads
-               if(indexOfPad){
-                 (nPlace[indexOfPad]<=AliTOFConstants::fgkNpadX) ? vol[4] = 1 : vol[4] = 2;
-                 (nPlace[indexOfPad]<=AliTOFConstants::fgkNpadX) ? vol[3] = nPlace[indexOfPad] : vol[3] = nPlace[indexOfPad] - AliTOFConstants::fgkNpadX;
-               }
-
-               // check if two sdigit are on the same pad; in that case we sum
-               // the two or more sdigits
-               if (hitMap->TestHit(vol) != kEmpty) {
-                 AliTOFSDigit *sdig = static_cast<AliTOFSDigit*>(hitMap->GetHit(vol));
-                 Int_t tdctime = (Int_t) digit[0];
-                 Int_t adccharge = (Int_t) digit[1];
-                 sdig->Update(tdctime,adccharge,tracknum);
-                 ntotalupdatesinEv++;
-                 ntotalupdates++;
-               } else {
-
-                 TOF->AddSDigit(tracknum, vol, digit);
 
+       // selection case for sdigitizing only hits in a given plate of a given sector
+       if(thereIsNotASelection || (vol[0]==fSelectedSector && vol[1]==fSelectedPlate)){
+         
+         vol[2] = tofHit->GetStrip();
+         vol[3] = tofHit->GetPadx();
+         vol[4] = tofHit->GetPadz();
+         
+         Bool_t dummy=((tracknum==previousTrack) && (vol[0]==previousSector) && (vol[1]==previousPlate) && (vol[2]==previousStrip));
+         
+         Bool_t isCloneOfThePrevious=dummy && ((vol[3]==previousPadX) && (vol[4]==previousPadZ));
+         
+         // much stronger check to be inserted here
+         
+         if(!isCloneOfThePrevious){
+           // update "previous" values
+           // in fact, we are yet in the future, so the present is past
+           previousTrack=tracknum;
+           previousSector=vol[0];
+           previousPlate=vol[1];
+           previousStrip=vol[2];
+           previousPadX=vol[3];
+           previousPadZ=vol[4];
+           
+           nselectedHits++;
+           nselectedHitsinEv++;
+           if (particle->GetFirstMother() < 0){
+             nHitsFromPrim++;
+           } // counts hits due to primary particles
+           
+           Float_t Xpad = tofHit->GetDx();
+           Float_t Zpad = tofHit->GetDz();
+           Float_t xStrip=AliTOFConstants::fgkXPad*(vol[3]-0.5-0.5*AliTOFConstants::fgkNpadX)+Xpad;
+           Float_t zStrip=AliTOFConstants::fgkZPad*(vol[4]-0.5-0.5*AliTOFConstants::fgkNpadZ)+Zpad;
+           Float_t geantTime = tofHit->GetTof(); // unit [s]
+           geantTime *= 1.e+09;  // conversion from [s] to [ns]
+           
+           //cout << "geantTime " << geantTime << " [ns]" << endl;
+           Int_t nActivatedPads = 0, nFiredPads = 0;
+           Bool_t isFired[4] = {kFALSE, kFALSE, kFALSE, kFALSE};
+           Float_t tofAfterSimul[4] = {0., 0., 0., 0.};
+           Float_t qInduced[4] = {0.,0.,0.,0.};
+           Int_t nPlace[4] = {0, 0, 0, 0};
+           Float_t averageTime = 0.;
+           SimulateDetectorResponse(zStrip,xStrip,geantTime,nActivatedPads,nFiredPads,isFired,nPlace,qInduced,tofAfterSimul,averageTime);
+           if(nFiredPads) {
+             for(Int_t indexOfPad=0; indexOfPad<nActivatedPads; indexOfPad++) {
+               if(isFired[indexOfPad]){ // the pad has fired
+                 Float_t timediff=geantTime-tofAfterSimul[indexOfPad];
+                 
+                 if(timediff>=0.2) nlargeTofDiff++;
+                 
+                 digit[0] = (Int_t) ((tofAfterSimul[indexOfPad]*1.e+03)/fTdcBin); // TDC bin number (each bin -> 25. ps)
+                 
+                 Float_t landauFactor = gRandom->Landau(fAdcMean, fAdcRms); 
+                 digit[1] = (Int_t) (qInduced[indexOfPad] * landauFactor); // ADC bins (each bin -> 0.25 (or 0.03) pC)
+                 
+                 // recalculate the volume only for neighbouring pads
                  if(indexOfPad){
-                   nnoisesdigits++;
-                   nnoisesdigitsinEv++;
-                 } else {
-                   nsignalsdigits++;
-                   nsignalsdigitsinEv++;
+                   (nPlace[indexOfPad]<=AliTOFConstants::fgkNpadX) ? vol[4] = 1 : vol[4] = 2;
+                   (nPlace[indexOfPad]<=AliTOFConstants::fgkNpadX) ? vol[3] = nPlace[indexOfPad] : vol[3] = nPlace[indexOfPad] - AliTOFConstants::fgkNpadX;
                  }
-                 ntotalsdigitsinEv++;  
-                 ntotalsdigits++;
-                 hitMap->SetHit(vol);
-               } // if (hitMap->TestHit(vol) != kEmpty)
-             } // if(isFired[indexOfPad])
-           } // end loop on nActivatedPads
-         } // if(nFiredPads) i.e. if some pads has fired
-       } // close if(!isCloneOfThePrevious)
+                 
+                 // check if two sdigit are on the same pad; in that case we sum
+                 // the two or more sdigits
+                 if (hitMap->TestHit(vol) != kEmpty) {
+                   AliTOFSDigit *sdig = static_cast<AliTOFSDigit*>(hitMap->GetHit(vol));
+                   Int_t tdctime = (Int_t) digit[0];
+                   Int_t adccharge = (Int_t) digit[1];
+                   sdig->Update(tdctime,adccharge,tracknum);
+                   ntotalupdatesinEv++;
+                   ntotalupdates++;
+                 } else {
+                   
+                   TOF->AddSDigit(tracknum, vol, digit);
+                   
+                   if(indexOfPad){
+                     nnoisesdigits++;
+                     nnoisesdigitsinEv++;
+                   } else {
+                     nsignalsdigits++;
+                     nsignalsdigitsinEv++;
+                   }
+                   ntotalsdigitsinEv++;  
+                   ntotalsdigits++;
+                   hitMap->SetHit(vol);
+                 } // if (hitMap->TestHit(vol) != kEmpty)
+               } // if(isFired[indexOfPad])
+             } // end loop on nActivatedPads
+           } // if(nFiredPads) i.e. if some pads has fired
+         } // close if(!isCloneOfThePrevious)
+       } // close the selection on sector and plate
       } // end loop on hits for the current track
     } // end loop on ntracks
-
+    
     delete hitMap;
       
     gAlice->TreeS()->Reset();
@@ -392,6 +414,22 @@ void AliTOFSDigitizer::Print(Option_t* opt)const
 
 }
 
+//__________________________________________________________________
+void AliTOFSDigitizer::SelectSectorAndPlate(Int_t sector, Int_t plate)
+{
+  Bool_t isaWrongSelection=(sector < 1) || (sector > AliTOFConstants::fgkNSectors) || (plate < 1) || (plate > AliTOFConstants::fgkNPlates);
+  if(isaWrongSelection){
+    cout << "You have selected an invalid value for sector or plate " << endl;
+    cout << "The correct range for sector is [1,"<< AliTOFConstants::fgkNSectors <<"]" << endl;
+    cout << "The correct range for plate  is [1,"<< AliTOFConstants::fgkNPlates  <<"]" << endl;
+    cout << "By default we continue sdigitizing all hits in all plates of all sectors" << endl;
+  } else {
+    fSelectedSector=sector;
+    fSelectedPlate =plate;
+    cout << "SDigitizing only hits in plate " << fSelectedPlate << " of the sector " << fSelectedSector << endl;
+  }
+}
+
 //__________________________________________________________________
 void AliTOFSDigitizer::SimulateDetectorResponse(Float_t z0, Float_t x0, Float_t geantTime, Int_t& nActivatedPads, Int_t& nFiredPads, Bool_t* isFired, Int_t* nPlace, Float_t* qInduced, Float_t* tofTime, Float_t& averageTime)
 {
index fd5491c..4eaf9a7 100644 (file)
@@ -37,6 +37,7 @@ public:
   Int_t GetFirstEvent()  const {return fEvent1;}
   Int_t GetSecondEvent() const {return fEvent2;}
   Int_t GetNEvents() const {return (fEvent2-fEvent1);}
+  void  SelectSectorAndPlate(Int_t sector, Int_t plate);
 
   // setters and getters for detector simulation
   // it summarizes all it is known about TOF strip 
@@ -107,6 +108,8 @@ private:
   Int_t   fEvent2;          // upper bound for events to sdigitize
   TF1     *ftail;           // pointer to formula for time with tail
   TString fHeadersFile;     // input file
+  Int_t fSelectedSector;    // sector number for sdigitization
+  Int_t fSelectedPlate ;    // plate  number for sdigitization
 
   // detector response simulation
   // Intrisic MRPC time resolution and pad (edge effect) parameters
@@ -152,7 +155,7 @@ private:
  protected:
 
 
-  ClassDef(AliTOFSDigitizer,1)  // creates TOF SDigits
+  ClassDef(AliTOFSDigitizer,2)  // creates TOF SDigits
 
 };
 
index c75cdac..8ab903a 100644 (file)
@@ -1,24 +1,36 @@
-Int_t AliTOFanalyzeSDigits(TString headersFile, Int_t iEvNum=0)
+Int_t AliTOFanalyzeSDigits(TString headersFile, Int_t ndump=15, Int_t iEvNum=0)
 {
   //
   // Analyzes the TOF sdigits and fills QA-histograms 
-  //
+  // report problems to pierella@bo.infn.it
   // iEvNum=0 means all events in the file
 
+  // Dynamically link some shared libs
+  if (gClassTable->GetID("AliRun") < 0) {
+    gROOT->LoadMacro("loadlibs.C");
+    loadlibs();
+  } else {
+    delete gAlice;
+    gAlice = 0;
+  }
+  
   Int_t rc=0;
 
-  TFile * file = (TFile*) gROOT->GetFile(headersFile.Data() ) ;
-
-  //File was not opened yet
-
-  if(file == 0){
-    if(headersFile.Contains("rfio"))
-      file =   TFile::Open(headersFile,"update") ;
-    else
-      file = new TFile(headersFile.Data(),"update") ;
-    gAlice = (AliRun *) file->Get("gAlice") ;
+  
+  TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(headersFile.Data());
+  if(file){
+    cout<<"headerFile already open \n";
   }
-
+  else {
+    if(!file)file=TFile::Open(headersFile.Data());
+  }
+  
+  // Get AliRun object from file
+  if (!gAlice) {
+    gAlice = (AliRun*)file->Get("gAlice");
+    if (gAlice) printf("AliRun object found on file\n");
+  }
+  
 
   if (iEvNum == 0) iEvNum = (Int_t) gAlice->TreeE()->GetEntries();
 
@@ -46,6 +58,7 @@ Int_t AliTOFanalyzeSDigits(TString headersFile, Int_t iEvNum=0)
   // ADC-TDC correlation
   TH2F *h2tdcVSadc = new TH2F("h2tdcVSadc","TDC [bin] VS ADC [bin]",500,0.,150000.,100,0.,3000.);
 
+  cout << "First " << ndump << " SDigits found in TOF TreeS branch have:" << endl;
 
   for (Int_t ievent = 0; ievent < iEvNum; ievent++) {
 
@@ -84,10 +97,16 @@ Int_t AliTOFanalyzeSDigits(TString headersFile, Int_t iEvNum=0)
       Bool_t isSDigitBad = (sector<1 || sector>18 || plate<1 || plate >5 || padz<1 || padz>2 || padx<1 || padx>48);
 
       if (isSDigitBad) {
-       cout << "<AliTOFanalyzeHits>  strange hit found" << endl;
+       cout << "<AliTOFanalyzeSDigits>  strange sdigit found" << endl;
        rc = 3;
        return rc;
       }
+      
+      if(k<ndump){
+       cout << k << "-th | " << "Sector " << sector << " | Plate " << plate << " | Strip " << strip << " | PadZ " << padz << " | PadX " << padx << endl;
+       cout << k << "-th | ADC " << firstADC << " [bin] | TDC " << firstTDC << " [bin]" << endl;
+       cout << "----------------------------------------------------"<< endl;
+      }
 
       // filling sdigit volume histos
       hsector->Fill(sector);
@@ -100,7 +119,7 @@ Int_t AliTOFanalyzeSDigits(TString headersFile, Int_t iEvNum=0)
       //cout << "firstTDC " << firstTDC << " firstADC " << firstADC << endl;
     }
   
-  }
+  } // end loop on events
 
   TFile *fout = new TFile("TOF_sdigitsQA.root","RECREATE");
   htdc->Write();
index 1ed178e..db2c6b4 100644 (file)
@@ -24,13 +24,21 @@ void AliTOFhits2sdigits(TString fileNameHits, Int_t firstEvent=0,Int_t nEvents=1
 
   // Create the TOF sdigitzer and sdigitize by default the first event
   // (in fact by default Int_t firstEvent=0,Int_t nEvents=1)
-  AliTOFSDigitizer *sdigitizer = new AliTOFSDigitizer(fileNameHits.Data(),fileNameHits.Data(),firstEvent,nEvents); // it is the same nevents numbering
+  AliTOFSDigitizer *sdigitizer = new AliTOFSDigitizer(fileNameHits.Data(),firstEvent,nEvents); // it is the same nevents numbering
   // scheme used by STEER/AliHits2SDigits.C
 
   // Activate this line if you want to print the parameters
   // used in sdigitization
   // sdigitizer->PrintParameters();
 
+  // e.g. Activate this line if you want to sdigitize only hits from plate 3
+  // in sector 15
+  // pay attention that sector must be in the range [1,18]
+  //                and plate  must be in the range [1,5]
+  // by default we sdigitize hits of all plates in all sectors
+  // sdigitizer->SelectSectorAndPlate(15,3);
+
+
   // performs sdigitization of the above events with "all" verbose option
   // "tim" option is also available for benchmarking only
   sdigitizer->Exec("all");