X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=EMCAL%2FAliEMCALSDigitizer.cxx;h=a483a5bae41d6f643ee9eaf666f02e01cb607e0c;hb=d7ea1ac71368bea2ccd3ca3d1d6c39c343f1b227;hp=24de52b9e2f76fb397744eaefece08c8016cdc68;hpb=4ce766eb4226c070c454ba72807fbf481119a5cb;p=u%2Fmrichter%2FAliRoot.git diff --git a/EMCAL/AliEMCALSDigitizer.cxx b/EMCAL/AliEMCALSDigitizer.cxx index 24de52b9e2f..a483a5bae41 100644 --- a/EMCAL/AliEMCALSDigitizer.cxx +++ b/EMCAL/AliEMCALSDigitizer.cxx @@ -16,7 +16,7 @@ /* $Id$ */ //_________________________________________________________________________ -// This is a TTask that makes SDigits out of Hits +// This is a Class that makes SDigits out of Hits // A Summable Digits is the sum of all hits originating // from one in one tower of the EMCAL // A threshold for assignment of the primary to SDigit is applied @@ -36,13 +36,13 @@ // User case: // root [0] AliEMCALSDigitizer * s = new AliEMCALSDigitizer("galice.root") // Warning in : object already instantiated -// root [1] s->ExecuteTask() +// root [1] s->Digitize() // // Makes SDigitis for all events stored in galice.root // root [2] s->SetPedestalParameter(0.001) // // One can change parameters of digitization // root [3] s->SetSDigitsBranch("Redestal 0.001") // // and write them into the new branch -// root [4] s->ExeciteTask("deb all tim") +// root [4] s->Digitize("deb all tim") // // available parameters: // deb - print # of produced SDigitis // deb all - print # and list of produced SDigits @@ -77,7 +77,7 @@ ClassImp(AliEMCALSDigitizer) //____________________________________________________________________________ AliEMCALSDigitizer::AliEMCALSDigitizer() - : TTask("",""), + : TNamed("",""), fA(0.),fB(0.),fECPrimThreshold(0.), fDefaultInit(kTRUE), fEventFolderName(0), @@ -85,7 +85,8 @@ AliEMCALSDigitizer::AliEMCALSDigitizer() fSDigitsInRun(0), fFirstEvent(0), fLastEvent(0), - fSampling(0.) + fSampling(0.), + fHits(0) { // ctor InitParameters(); @@ -94,7 +95,7 @@ AliEMCALSDigitizer::AliEMCALSDigitizer() //____________________________________________________________________________ AliEMCALSDigitizer::AliEMCALSDigitizer(const char * alirunFileName, const char * eventFolderName) - : TTask("EMCAL"+AliConfig::Instance()->GetSDigitizerTaskName(), alirunFileName), + : TNamed("EMCALSDigitizer", alirunFileName), fA(0.),fB(0.),fECPrimThreshold(0.), fDefaultInit(kFALSE), fEventFolderName(eventFolderName), @@ -102,7 +103,8 @@ AliEMCALSDigitizer::AliEMCALSDigitizer(const char * alirunFileName, fSDigitsInRun(0), fFirstEvent(0), fLastEvent(0), - fSampling(0.) + fSampling(0.), + fHits(0) { // ctor Init(); @@ -113,7 +115,7 @@ AliEMCALSDigitizer::AliEMCALSDigitizer(const char * alirunFileName, //____________________________________________________________________________ AliEMCALSDigitizer::AliEMCALSDigitizer(const AliEMCALSDigitizer & sd) - : TTask(sd.GetName(),sd.GetTitle()), + : TNamed(sd.GetName(),sd.GetTitle()), fA(sd.fA), fB(sd.fB), fECPrimThreshold(sd.fECPrimThreshold), @@ -123,24 +125,34 @@ AliEMCALSDigitizer::AliEMCALSDigitizer(const AliEMCALSDigitizer & sd) fSDigitsInRun(sd.fSDigitsInRun), fFirstEvent(sd.fFirstEvent), fLastEvent(sd.fLastEvent), - fSampling(sd.fSampling) + fSampling(sd.fSampling), + fHits(0) { //cpy ctor } +//_____________________________________________________________________ +AliEMCALSDigitizer& AliEMCALSDigitizer::operator = (const AliEMCALSDigitizer& source) +{ // assignment operator; use copy ctor + if (&source == this) return *this; + + new (this) AliEMCALSDigitizer(source); + return *this; +} //____________________________________________________________________________ AliEMCALSDigitizer::~AliEMCALSDigitizer() { //dtor - AliLoader *emcalLoader=0; - if ((emcalLoader = AliRunLoader::Instance()->GetDetectorLoader("EMCAL"))) - emcalLoader->CleanSDigitizer(); + + if(fHits){ + fHits->Clear(); + delete fHits; + } } //____________________________________________________________________________ void AliEMCALSDigitizer::Init(){ - // Initialization: open root-file, allocate arrays for hits and sdigits, - // attach task SDigitizer to the list of EMCAL tasks + // Initialization: open root-file, allocate arrays for hits and sdigits // // Initialization can not be done in the default constructor //============================================================= YS @@ -155,9 +167,6 @@ void AliEMCALSDigitizer::Init(){ return ; } - emcalLoader->PostSDigitizer(this); - emcalLoader->GetSDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE); - } //____________________________________________________________________________ @@ -212,148 +221,178 @@ void AliEMCALSDigitizer::InitParameters() } //____________________________________________________________________________ -void AliEMCALSDigitizer::Exec(Option_t *option) +void AliEMCALSDigitizer::Digitize(Option_t *option) { - // Collects all hit of the same tower into digits - TString o(option); o.ToUpper(); - if (strstr(option, "print") ) { - - AliDebug(2,Form("Print: \n------------------- %s -------------\n",GetName())); - AliDebug(2,Form(" fInit %i\n", int(fInit))); - AliDebug(2,Form(" fFirstEvent %i\n", fFirstEvent)); - AliDebug(2,Form(" fLastEvent %i\n", fLastEvent)); - AliDebug(2,Form(" Writing SDigits to branch with title %s\n", fEventFolderName.Data())); - AliDebug(2,Form(" with digitization parameters A = %f\n", fA)); - AliDebug(2,Form(" B = %f\n", fB)); - AliDebug(2,Form(" Threshold for EC Primary assignment = %f\n", fECPrimThreshold)); - AliDebug(2,Form(" Sampling = %f\n", fSampling)); - AliDebug(2,Form("---------------------------------------------------\n")); - - return ; - } - - - if(strstr(option,"tim")) - gBenchmark->Start("EMCALSDigitizer"); - - AliRunLoader *rl = AliRunLoader::Instance(); - AliEMCALLoader *emcalLoader = dynamic_cast(rl->GetDetectorLoader("EMCAL")); - - //switch off reloading of this task while getting event - if (!fInit) { // to prevent overwrite existing file - AliError( Form("Give a version name different from %s", fEventFolderName.Data()) ) ; - return ; - } - - if (fLastEvent == -1) - fLastEvent = rl->GetNumberOfEvents() - 1 ; - else { - fLastEvent = TMath::Min(fLastEvent, rl->GetNumberOfEvents()-1); - } - Int_t nEvents = fLastEvent - fFirstEvent + 1; - - Int_t ievent; - Float_t energy=0.; // de * fSampling - 23-nov-04 - rl->LoadKinematics(); - rl->LoadHits("EMCAL"); - - for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) { - rl->GetEvent(ievent); - TTree * treeS = emcalLoader->TreeS(); - if ( !treeS ) { - emcalLoader->MakeSDigitsContainer(); - treeS = emcalLoader->TreeS(); - } - TClonesArray * hits = emcalLoader->Hits() ; - TClonesArray * sdigits = emcalLoader->SDigits() ; - sdigits->Clear(); - - Int_t nSdigits = 0 ; - Int_t i; - - AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance(); - for ( i = 0 ; i < hits->GetEntries() ; i++ ) { - AliEMCALHit * hit = dynamic_cast(hits->At(i)) ; - AliEMCALDigit * curSDigit = 0 ; - AliEMCALDigit * sdigit = 0 ; - Bool_t newsdigit = kTRUE; - - // hit->GetId() - Absolute Id number EMCAL segment - if(geom->CheckAbsCellId(hit->GetId())) { // was IsInECA(hit->GetId()) - energy = hit->GetEnergy() * fSampling; // 23-nov-04 - if(energy > fECPrimThreshold ) - // Assign primary number only if deposited energy is significant - curSDigit = new AliEMCALDigit( hit->GetPrimary(), - hit->GetIparent(), hit->GetId(), - Digitize(energy), hit->GetTime(),kFALSE, - -1, energy ) ; - else - curSDigit = new AliEMCALDigit( -1, - -1, - hit->GetId(), - Digitize(energy), hit->GetTime(),kFALSE, - -1, energy ) ; - } else { - Warning("Exec"," abs id %i is bad \n", hit->GetId()); - newsdigit = kFALSE; - curSDigit = 0; - } - - if(curSDigit != 0){ - for(Int_t check= 0; check < nSdigits ; check++) { - sdigit = dynamic_cast(sdigits->At(check)) ; - - if( sdigit->GetId() == curSDigit->GetId()) { // Are we in the same ECAL tower ? - *sdigit = *sdigit + *curSDigit; - newsdigit = kFALSE; - } + // Collects all hit of the same tower into digits + TString o(option); o.ToUpper(); + if (strstr(option, "print") ) { + + AliDebug(2,Form("Print: \n------------------- %s -------------\n",GetName())); + AliDebug(2,Form(" fInit %i\n", int(fInit))); + AliDebug(2,Form(" fFirstEvent %i\n", fFirstEvent)); + AliDebug(2,Form(" fLastEvent %i\n", fLastEvent)); + AliDebug(2,Form(" Writing SDigits to branch with title %s\n", fEventFolderName.Data())); + AliDebug(2,Form(" with digitization parameters A = %f\n", fA)); + AliDebug(2,Form(" B = %f\n", fB)); + AliDebug(2,Form(" Threshold for EC Primary assignment = %f\n", fECPrimThreshold)); + AliDebug(2,Form(" Sampling = %f\n", fSampling)); + AliDebug(2,Form("---------------------------------------------------\n")); + + return ; } - } - if (newsdigit) { - new((*sdigits)[nSdigits]) AliEMCALDigit(*curSDigit); - nSdigits++ ; - } - delete curSDigit ; - } // loop over all hits (hit = deposited energy/entering particle) - sdigits->Sort() ; - - nSdigits = sdigits->GetEntriesFast() ; - fSDigitsInRun += nSdigits ; - - for (i = 0 ; i < sdigits->GetEntriesFast() ; i++) { - AliEMCALDigit * sdigit = dynamic_cast(sdigits->At(i)) ; - sdigit->SetIndexInList(i) ; - } - - // Now write SDigits - - Int_t bufferSize = 32000 ; - TBranch * sdigitsBranch = treeS->GetBranch("EMCAL"); - if (sdigitsBranch) - sdigitsBranch->SetAddress(&sdigits); - else - treeS->Branch("EMCAL",&sdigits,bufferSize); - - treeS->Fill(); - - emcalLoader->WriteSDigits("OVERWRITE"); - - //NEXT - SDigitizer - emcalLoader->WriteSDigitizer("OVERWRITE"); // why in event cycle ? - - if(strstr(option,"deb")) - PrintSDigits(option) ; - } - - Unload(); - - emcalLoader->GetSDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE); + + + if(strstr(option,"tim")) + gBenchmark->Start("EMCALSDigitizer"); + + AliRunLoader *rl = AliRunLoader::Instance(); + AliEMCALLoader *emcalLoader = dynamic_cast(rl->GetDetectorLoader("EMCAL")); + + if (!fInit) { // to prevent overwrite existing file + AliError( Form("Give a version name different from %s", fEventFolderName.Data()) ) ; + return ; + } + + if (fLastEvent == -1) + fLastEvent = rl->GetNumberOfEvents() - 1 ; + else { + fLastEvent = TMath::Min(fLastEvent, rl->GetNumberOfEvents()-1); + } + Int_t nEvents = fLastEvent - fFirstEvent + 1; + + Float_t energy=0.; // de * fSampling - 23-nov-04 + rl->LoadKinematics(); + rl->LoadHits("EMCAL"); - if(strstr(option,"tim")){ - gBenchmark->Stop("EMCALSDigitizer"); - printf("\n Exec: took %f seconds for SDigitizing %f seconds per event\n", - gBenchmark->GetCpuTime("EMCALSDigitizer"), gBenchmark->GetCpuTime("EMCALSDigitizer")/nEvents ) ; - } + Int_t ievent; + for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) { + rl->GetEvent(ievent); + TTree * treeS = emcalLoader->TreeS(); + if ( !treeS ) { + emcalLoader->MakeSDigitsContainer(); + treeS = emcalLoader->TreeS(); + } + + TClonesArray * sdigits = emcalLoader->SDigits() ; + sdigits->Clear("C"); + + Int_t nSdigits = 0 ; + Int_t iHit, iTrack, iSDigit; + + AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance(); + + TTree *treeH = emcalLoader->TreeH(); + if (treeH) { + Int_t nTrack = treeH->GetEntries(); // TreeH has array of hits for every primary + TBranch * branchH = treeH->GetBranch("EMCAL"); + //if(fHits)fHits->Clear(); + branchH->SetAddress(&fHits); + for (iTrack = 0; iTrack < nTrack; iTrack++) { + branchH->GetEntry(iTrack); + + if(fHits){ + + Int_t nHit = fHits->GetEntriesFast(); + for(iHit = 0; iHit< nHit;iHit++){ + + AliEMCALHit * hit = dynamic_cast(fHits->At(iHit)) ; + AliEMCALDigit * curSDigit = 0 ; + AliEMCALDigit * sdigit = 0 ; + Bool_t newsdigit = kTRUE; + + // hit->GetId() - Absolute Id number EMCAL segment + if(hit){ + if(geom->CheckAbsCellId(hit->GetId())) { // was IsInECA(hit->GetId()) + energy = hit->GetEnergy() * fSampling; // 23-nov-04 + if(energy > fECPrimThreshold ) + // Assign primary number only if deposited energy is significant + curSDigit = new AliEMCALDigit(hit->GetPrimary(), + hit->GetIparent(), hit->GetId(), + Digitize(energy), hit->GetTime(),kFALSE, + -1, 0,0,energy ) ; + else + curSDigit = new AliEMCALDigit(-1, + -1, + hit->GetId(), + Digitize(energy), hit->GetTime(),kFALSE, + -1, 0,0,energy ) ; + } else { + Warning("Digitize"," abs id %i is bad \n", hit->GetId()); + newsdigit = kFALSE; + curSDigit = 0; + } + + if(curSDigit != 0){ + for(Int_t check= 0; check < nSdigits ; check++) { + sdigit = dynamic_cast(sdigits->At(check)) ; + if(sdigit){ + if( sdigit->GetId() == curSDigit->GetId()) { // Are we in the same ECAL tower ? + *sdigit = *sdigit + *curSDigit; + newsdigit = kFALSE; + } + }// sdigit exists + else { + AliWarning("Sdigit do not exist"); + newsdigit = kFALSE; + }// sdigit does not exist + }//sdigit loop + }// currsdigit exists + + if (newsdigit) { + new((*sdigits)[nSdigits]) AliEMCALDigit(*curSDigit); + nSdigits++ ; + } + delete curSDigit ; + + }// hit exists + else AliFatal("Hit is NULL!"); + + } // loop over all hits (hit = deposited energy/entering particle) + + }//fHits is not NULL + else AliFatal("fHit is NULL!"); + + sdigits->Sort() ; + + nSdigits = sdigits->GetEntriesFast() ; + fSDigitsInRun += nSdigits ; + + for (iSDigit = 0 ; iSDigit < sdigits->GetEntriesFast() ; iSDigit++) { + AliEMCALDigit * sdigit = dynamic_cast(sdigits->At(iSDigit)) ; + if(sdigit)sdigit->SetIndexInList(iSDigit) ; + else AliFatal("sdigit is NULL!"); + } + if(fHits)fHits->Clear(); + }//track loop + }// tree exists + + // Now write SDigits + + Int_t bufferSize = 32000 ; + TBranch * sdigitsBranch = treeS->GetBranch("EMCAL"); + if (sdigitsBranch) + sdigitsBranch->SetAddress(&sdigits); + else + treeS->Branch("EMCAL",&sdigits,bufferSize); + + treeS->Fill(); + + emcalLoader->WriteSDigits("OVERWRITE"); + + //NEXT - SDigitizer + //emcalLoader->WriteSDigitizer("OVERWRITE"); // why in event cycle ? + + if(strstr(option,"deb")) + PrintSDigits(option) ; + } + + Unload(); + + if(strstr(option,"tim")){ + gBenchmark->Stop("EMCALSDigitizer"); + printf("\n Digitize: took %f seconds for SDigitizing %f seconds per event\n", + gBenchmark->GetCpuTime("EMCALSDigitizer"), gBenchmark->GetCpuTime("EMCALSDigitizer")/nEvents ) ; + } } //__________________________________________________________________ @@ -432,42 +471,48 @@ Bool_t AliEMCALSDigitizer::operator==( AliEMCALSDigitizer const &sd )const else return kFALSE ; } - //__________________________________________________________________ void AliEMCALSDigitizer::PrintSDigits(Option_t * option) { //Prints list of digits produced at the current pass of AliEMCALDigitizer - - AliEMCALLoader *rl = dynamic_cast(AliRunLoader::Instance()->GetDetectorLoader("EMCAL")); - const TClonesArray * sdigits = rl->SDigits() ; - printf("\n") ; - printf("event %i", rl->GetRunLoader()->GetEventNumber()); - printf(" Number of entries in SDigits list %i", sdigits->GetEntriesFast()); - if(strstr(option,"all")||strstr(option,"EMC")){ + AliEMCALLoader *rl = dynamic_cast(AliRunLoader::Instance()->GetDetectorLoader("EMCAL")); + if(rl){ + const TClonesArray * sdigits = rl->SDigits() ; - //loop over digits - AliEMCALDigit * digit; - printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ; - Int_t index = 0; - Float_t isum = 0.; - char * tempo = new char[8192]; - for (index = 0 ; index < sdigits->GetEntries() ; index++) { - digit = dynamic_cast( sdigits->At(index) ) ; - sprintf(tempo, "\n%6d %8f %6.5e %4d %2d :", - digit->GetId(), digit->GetAmplitude(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; - printf("%s",tempo); - isum += digit->GetAmplitude(); + printf("\n") ; + printf("event %i", rl->GetRunLoader()->GetEventNumber()); + printf(" Number of entries in SDigits list %i", sdigits->GetEntriesFast()); + if(strstr(option,"all")||strstr(option,"EMC")){ - Int_t iprimary; - for (iprimary=0; iprimaryGetNprimary(); iprimary++) { - sprintf(tempo, "%d ",digit->GetPrimary(iprimary+1) ) ; - printf("%s",tempo); - } - } - delete [] tempo ; - printf("\n** Sum %2.3f : %10.3f GeV/c **\n ", isum, Calibrate(isum)); - } else printf("\n"); + //loop over digits + AliEMCALDigit * digit; + printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ; + Int_t index = 0; + Float_t isum = 0.; + const Int_t bufferSize = 8192; + char * tempo = new char[bufferSize]; + for (index = 0 ; index < sdigits->GetEntries() ; index++) { + digit = dynamic_cast( sdigits->At(index) ) ; + if(digit){ + snprintf(tempo, bufferSize,"\n%6d %8f %6.5e %4d %2d :", + digit->GetId(), digit->GetAmplitude(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; + printf("%s",tempo); + isum += digit->GetAmplitude(); + + Int_t iprimary; + for (iprimary=0; iprimaryGetNprimary(); iprimary++) { + snprintf(tempo,bufferSize, "%d ",digit->GetPrimary(iprimary+1) ) ; + printf("%s",tempo); + } + } //sdigit exists + else AliFatal("SDigit is NULL!"); + }//loop + delete [] tempo ; + printf("\n** Sum %2.3f : %10.3f GeV/c **\n ", isum, Calibrate(isum)); + } else printf("\n"); + } + else AliFatal("EMCALLoader is NULL!"); } //____________________________________________________________________________ @@ -475,12 +520,15 @@ void AliEMCALSDigitizer::Unload() const { // Unload Hits and SDigits from the folder AliEMCALLoader *rl = dynamic_cast(AliRunLoader::Instance()->GetDetectorLoader("EMCAL")); - rl->UnloadHits() ; - rl->UnloadSDigits() ; + if(rl){ + rl->UnloadHits() ; + rl->UnloadSDigits() ; + } + else AliFatal("EMCALLoader is NULL!"); } //____________________________________________________________________________ void AliEMCALSDigitizer::Browse(TBrowser* b) { - TTask::Browse(b); + TNamed::Browse(b); }