X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PHOS%2FAliPHOSClusterizerv1.cxx;h=968c6c425a67a8e1a01c47dfbbd6f81a4a09922c;hb=28871337c114f0a762da76ff09de652633e3dde1;hp=d001ec36f7920d53da528e322ca96f407699783f;hpb=fb43ada4db956ceb59baa593557bff7082e4bfc2;p=u%2Fmrichter%2FAliRoot.git diff --git a/PHOS/AliPHOSClusterizerv1.cxx b/PHOS/AliPHOSClusterizerv1.cxx index d001ec36f79..968c6c425a6 100644 --- a/PHOS/AliPHOSClusterizerv1.cxx +++ b/PHOS/AliPHOSClusterizerv1.cxx @@ -15,6 +15,44 @@ /* $Id$ */ +/* History of cvs commits: + * + * $Log$ + * Revision 1.89 2006/03/13 14:05:42 kharlov + * Calibration objects for EMC and CPV + * + * Revision 1.88 2006/01/11 08:54:52 hristov + * Additional protection in case no calibration entry was found + * + * Revision 1.87 2005/11/22 08:46:43 kharlov + * Updated to new CDB (Boris Polichtchouk) + * + * Revision 1.86 2005/11/14 21:52:43 hristov + * Coding conventions + * + * Revision 1.85 2005/09/27 16:08:08 hristov + * New version of CDB storage framework (A.Colla) + * + * Revision 1.84 2005/09/21 10:02:47 kharlov + * Reading calibration from CDB (Boris Polichtchouk) + * + * Revision 1.82 2005/09/02 15:43:13 kharlov + * Add comments in GetCalibrationParameters and Calibrate + * + * Revision 1.81 2005/09/02 14:32:07 kharlov + * Calibration of raw data + * + * Revision 1.80 2005/08/24 15:31:36 kharlov + * Setting raw digits flag + * + * Revision 1.79 2005/07/25 15:53:53 kharlov + * Read raw data + * + * Revision 1.78 2005/05/28 14:19:04 schutz + * Compilation warnings fixed by T.P. + * + */ + //*-- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (SUBATECH & Kurchatov Institute) ////////////////////////////////////////////////////////////////////////////// // Clusterization class. Performs clusterization (collects neighbouring active cells) and @@ -53,6 +91,7 @@ // --- Standard library --- // --- AliRoot header files --- +#include "AliLog.h" #include "AliPHOSGetter.h" #include "AliPHOSGeometry.h" #include "AliPHOSClusterizerv1.h" @@ -60,6 +99,10 @@ #include "AliPHOSCpvRecPoint.h" #include "AliPHOSDigit.h" #include "AliPHOSDigitizer.h" +#include "AliPHOSCalibrationDB.h" +#include "AliCDBManager.h" +#include "AliCDBStorage.h" +#include "AliCDBEntry.h" ClassImp(AliPHOSClusterizerv1) @@ -88,11 +131,7 @@ AliPHOSClusterizerv1::AliPHOSClusterizerv1(const TString alirunFileName, const T { // dtor - // delete fPedestals ; -// delete fGains ; } - - //____________________________________________________________________________ const TString AliPHOSClusterizerv1::BranchName() const { @@ -100,17 +139,36 @@ const TString AliPHOSClusterizerv1::BranchName() const } //____________________________________________________________________________ -Float_t AliPHOSClusterizerv1::Calibrate(Int_t amp, Int_t absId) const +Float_t AliPHOSClusterizerv1::Calibrate(Int_t amp, Int_t absId) { - //To be replaced later by the method, reading individual parameters from the database - // if(fCalibrationDB) - // return fCalibrationDB->Calibrate(amp,absId) ; - // else{ //simulation - if(absId <= fEmcCrystals) //calibrate as EMC - return fADCpedestalEmc + amp*fADCchanelEmc ; - else //calibrate as CPV - return fADCpedestalCpv+ amp*fADCchanelCpv ; - // } + // Convert digitized amplitude into energy. + // Calibration parameters are taken from calibration data base for raw data, + // or from digitizer parameters for simulated data. + + if(fCalibData){ + Int_t relId[4]; + AliPHOSGetter *gime = AliPHOSGetter::Instance(); + gime->PHOSGeometry()->AbsToRelNumbering(absId,relId) ; + Int_t module = relId[0]; + Int_t column = relId[3]; + Int_t row = relId[2]; + if(absId <= fEmcCrystals) { //calibrate as EMC + fADCchanelEmc = fCalibData->GetADCchannelEmc (module,column,row); + fADCpedestalEmc = fCalibData->GetADCpedestalEmc(module,column,row); + return fADCpedestalEmc + amp*fADCchanelEmc ; + } + else { //calibrate as CPV + fADCchanelCpv = fCalibData->GetADCchannelCpv (module,column,row); + fADCpedestalCpv = fCalibData->GetADCpedestalCpv(module,column,row); + return fADCpedestalCpv + amp*fADCchanelCpv ; + } + } + else{ //simulation + if(absId <= fEmcCrystals) //calibrate as EMC + return fADCpedestalEmc + amp*fADCchanelEmc ; + else //calibrate as CPV + return fADCpedestalCpv+ amp*fADCchanelCpv ; + } } //____________________________________________________________________________ @@ -129,20 +187,27 @@ void AliPHOSClusterizerv1::Exec(Option_t *option) return ; } - AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ; + GetCalibrationParameters() ; + + AliPHOSGetter * gime = AliPHOSGetter::Instance() ; + if (fRawReader == 0) + gime->SetRawDigits(kFALSE); + else + gime->SetRawDigits(kTRUE); if (fLastEvent == -1) fLastEvent = gime->MaxEvent() - 1 ; else fLastEvent = TMath::Min(fFirstEvent, gime->MaxEvent()); // one event at the time Int_t nEvents = fLastEvent - fFirstEvent + 1; - + Int_t ievent ; for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) { - gime->Event(ievent, "D"); - - GetCalibrationParameters() ; - + if (fRawReader == 0) + gime->Event(ievent ,"D"); // Read digits from simulated data + else + gime->Event(fRawReader,"W"); // Read digits from raw data + fNumberOfEmcClusters = fNumberOfCpvClusters = 0 ; MakeClusters() ; @@ -160,13 +225,14 @@ void AliPHOSClusterizerv1::Exec(Option_t *option) fRecPointsInRun += gime->CpvRecPoints()->GetEntriesFast() ; } - Unload(); + if(fWrite) //do not unload in "on flight" mode + Unload(); if(strstr(option,"tim")){ gBenchmark->Stop("PHOSClusterizer"); - Info("Exec", " took %f seconds for Clusterizing %f seconds per event \n", + AliInfo(Form("took %f seconds for Clusterizing %f seconds per event \n", gBenchmark->GetCpuTime("PHOSClusterizer"), - gBenchmark->GetCpuTime("PHOSClusterizer")/nEvents ) ; + gBenchmark->GetCpuTime("PHOSClusterizer")/nEvents )) ; } } @@ -267,20 +333,29 @@ Bool_t AliPHOSClusterizerv1::FindFit(AliPHOSEmcRecPoint * emcRP, AliPHOSDigit ** //____________________________________________________________________________ void AliPHOSClusterizerv1::GetCalibrationParameters() { + // Set calibration parameters: + // if calibration database exists, they are read from database, + // otherwise, they are taken from digitizer. + // + // It is a user responsilibity to open CDB before reconstruction, for example: + // AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB"); AliPHOSGetter * gime = AliPHOSGetter::Instance(); - - if ( !gime->Digitizer() ) - gime->LoadDigitizer(); - AliPHOSDigitizer * dig = gime->Digitizer(); - - fADCchanelEmc = dig->GetEMCchannel() ; - fADCpedestalEmc = dig->GetEMCpedestal(); - - fADCchanelCpv = dig->GetCPVchannel() ; - fADCpedestalCpv = dig->GetCPVpedestal() ; + // fCalibData = new AliPHOSCalibData(gAlice->GetRunNumber()); //original + + fCalibData = new AliPHOSCalibData(-1); //use AliCDBManager's run number - // fCalibrationDB = gime->CalibrationDB(); + if(!fCalibData) + { + if ( !gime->Digitizer() ) + gime->LoadDigitizer(); + AliPHOSDigitizer * dig = gime->Digitizer(); + fADCchanelEmc = dig->GetEMCchannel() ; + fADCpedestalEmc = dig->GetEMCpedestal(); + + fADCchanelCpv = dig->GetCPVchannel() ; + fADCpedestalCpv = dig->GetCPVpedestal() ; + } } //____________________________________________________________________________ @@ -288,18 +363,21 @@ void AliPHOSClusterizerv1::Init() { // Make all memory allocations which can not be done in default constructor. // Attach the Clusterizer task to the list of PHOS tasks - - AliPHOSGetter* gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data()); + + AliPHOSGetter* gime = AliPHOSGetter::Instance() ; + if(!gime) + gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data()); AliPHOSGeometry * geom = gime->PHOSGeometry(); - + fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ; if(!gMinuit) gMinuit = new TMinuit(100); - if ( !gime->Clusterizer() ) + if ( !gime->Clusterizer() ) { gime->PostClusterizer(this); + } } //____________________________________________________________________________ @@ -314,7 +392,10 @@ void AliPHOSClusterizerv1::InitParameters() fEmcLocMaxCut = 0.03 ; fCpvLocMaxCut = 0.03 ; - + + fEmcMinE = 0.01 ; + fCpvMinE = 0.0 ; + fW0 = 4.5 ; fW0CPV = 4.0 ; @@ -324,6 +405,10 @@ void AliPHOSClusterizerv1::InitParameters() fRecPointsInRun = 0 ; + fWrite = kTRUE ; + + fCalibData = 0 ; + SetEventRange(0,-1) ; } @@ -370,8 +455,21 @@ Int_t AliPHOSClusterizerv1::AreNeighbours(AliPHOSDigit * d1, AliPHOSDigit * d2)c return rv ; } +//____________________________________________________________________________ +void AliPHOSClusterizerv1::CleanDigits(TClonesArray * digits){ + for(Int_t i=0; iGetEntriesFast(); i++){ + AliPHOSDigit * digit = static_cast(digits->At(i)) ; + Float_t cut = IsInEmc(digit) ? fEmcMinE : fCpvMinE ; + if(Calibrate(digit->GetAmp(),digit->GetId()) < cut) + digits->RemoveAt(i) ; + } + digits->Compress() ; + for (Int_t i = 0 ; i < digits->GetEntriesFast() ; i++) { + AliPHOSDigit *digit = static_cast( digits->At(i) ) ; + digit->SetIndexInList(i) ; + } - +} //____________________________________________________________________________ Bool_t AliPHOSClusterizerv1::IsInEmc(AliPHOSDigit * digit) const { @@ -423,47 +521,60 @@ void AliPHOSClusterizerv1::WriteRecPoints() TObjArray * emcRecPoints = gime->EmcRecPoints() ; TObjArray * cpvRecPoints = gime->CpvRecPoints() ; TClonesArray * digits = gime->Digits() ; - - TTree * treeR = gime->TreeR(); - + + Int_t index ; - //Evaluate position, dispersion and other RecPoint properties... - for(index = 0; index < emcRecPoints->GetEntries(); index++) - dynamic_cast( emcRecPoints->At(index) )->EvalAll(fW0,digits) ; - + //Evaluate position, dispersion and other RecPoint properties.. + Int_t nEmc = emcRecPoints->GetEntriesFast(); + for(index = 0; index < nEmc; index++){ + AliPHOSEmcRecPoint * rp = dynamic_cast( emcRecPoints->At(index) ); + rp->Purify(fEmcMinE) ; + if(rp->GetMultiplicity()>0.) //If this RP is not empty + rp->EvalAll(fW0,digits) ; + else{ + emcRecPoints->RemoveAt(index) ; + delete rp ; + } + } + emcRecPoints->Compress() ; emcRecPoints->Sort() ; - for(index = 0; index < emcRecPoints->GetEntries(); index++) + // emcRecPoints->Expand(emcRecPoints->GetEntriesFast()) ; + for(index = 0; index < emcRecPoints->GetEntries(); index++){ dynamic_cast( emcRecPoints->At(index) )->SetIndexInList(index) ; - - emcRecPoints->Expand(emcRecPoints->GetEntriesFast()) ; + } //Now the same for CPV - for(index = 0; index < cpvRecPoints->GetEntries(); index++) - dynamic_cast( cpvRecPoints->At(index) )->EvalAll(digits) ; - + for(index = 0; index < cpvRecPoints->GetEntries(); index++){ + AliPHOSCpvRecPoint * rp = dynamic_cast( cpvRecPoints->At(index) ); + rp->EvalAll(fW0CPV,digits) ; + } cpvRecPoints->Sort() ; for(index = 0; index < cpvRecPoints->GetEntries(); index++) - dynamic_cast( cpvRecPoints->At(index) )->SetIndexInList(index) ; + dynamic_cast( cpvRecPoints->At(index) )->SetIndexInList(index) ; cpvRecPoints->Expand(cpvRecPoints->GetEntriesFast()) ; - Int_t bufferSize = 32000 ; - Int_t splitlevel = 0 ; - - //First EMC - TBranch * emcBranch = treeR->Branch("PHOSEmcRP","TObjArray",&emcRecPoints,bufferSize,splitlevel); - emcBranch->SetTitle(BranchName()); + if(fWrite){ //We write TreeR + TTree * treeR = gime->TreeR(); - //Now CPV branch - TBranch * cpvBranch = treeR->Branch("PHOSCpvRP","TObjArray",&cpvRecPoints,bufferSize,splitlevel); - cpvBranch->SetTitle(BranchName()); + Int_t bufferSize = 32000 ; + Int_t splitlevel = 0 ; - emcBranch ->Fill() ; - cpvBranch ->Fill() ; - - gime->WriteRecPoints("OVERWRITE"); - gime->WriteClusterizer("OVERWRITE"); + //First EMC + TBranch * emcBranch = treeR->Branch("PHOSEmcRP","TObjArray",&emcRecPoints,bufferSize,splitlevel); + emcBranch->SetTitle(BranchName()); + + //Now CPV branch + TBranch * cpvBranch = treeR->Branch("PHOSCpvRP","TObjArray",&cpvRecPoints,bufferSize,splitlevel); + cpvBranch->SetTitle(BranchName()); + + emcBranch ->Fill() ; + cpvBranch ->Fill() ; + + gime->WriteRecPoints("OVERWRITE"); + gime->WriteClusterizer("OVERWRITE"); + } } //____________________________________________________________________________ @@ -483,6 +594,10 @@ void AliPHOSClusterizerv1::MakeClusters() TClonesArray * digits = gime->Digits() ; + //Remove digits below threshold + CleanDigits(digits) ; + + TClonesArray * digitsC = static_cast( digits->Clone() ) ; @@ -625,6 +740,9 @@ void AliPHOSClusterizerv1::MakeUnfolding() fNumberOfEmcClusters -- ; numberofNotUnfolded-- ; } + else{ + emcRecPoint->SetNExMax(1) ; //Only one local maximum + } delete[] maxAt ; delete[] maxAtEnergy ; @@ -705,6 +823,7 @@ void AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc, Bool_t rv = FindFit(iniEmc, maxAt, maxAtEnergy, nPar, fitparameters) ; if( !rv ) { // Fit failed, return and remove cluster + iniEmc->SetNExMax(-1) ; delete[] fitparameters ; return ; } @@ -766,6 +885,7 @@ void AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc, (*emcRecPoints)[fNumberOfEmcClusters] = new AliPHOSEmcRecPoint("") ; emcRP = dynamic_cast( emcRecPoints->At(fNumberOfEmcClusters) ) ; fNumberOfEmcClusters++ ; + emcRP->SetNExMax((Int_t)nPar/3) ; } else{//create new entries in fCpvRecPoints if(fNumberOfCpvClusters >= cpvRecPoints->GetSize()) @@ -893,7 +1013,7 @@ void AliPHOSClusterizerv1::UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Dou } //____________________________________________________________________________ -void AliPHOSClusterizerv1::Print()const +void AliPHOSClusterizerv1::Print(const Option_t *)const { // Print clusterizer parameters @@ -922,7 +1042,7 @@ void AliPHOSClusterizerv1::Print()const else message = " AliPHOSClusterizerv1 not initialized " ; - Info("Print", message.Data(), + AliInfo(Form("%s, %s %s %s %s %s %s %s %s %s %s", message.Data(), taskName.Data(), GetTitle(), taskName.Data(), @@ -932,7 +1052,7 @@ void AliPHOSClusterizerv1::Print()const fW0, fCpvClusteringThreshold, fCpvLocMaxCut, - fW0CPV ) ; + fW0CPV )) ; } @@ -946,9 +1066,10 @@ void AliPHOSClusterizerv1::PrintRecPoints(Option_t * option) TObjArray * emcRecPoints = gime->EmcRecPoints() ; TObjArray * cpvRecPoints = gime->CpvRecPoints() ; - printf("\nevent %d \n", gAlice->GetEvNumber()) ; - printf(" Found %d EMC RecPoints and %d CPV RecPoints \n", - emcRecPoints->GetEntriesFast(),cpvRecPoints->GetEntriesFast()) ; + AliInfo(Form("\nevent %d \n Found %d EMC RecPoints and %d CPV RecPoints", + gAlice->GetEvNumber(), + emcRecPoints->GetEntriesFast(), + cpvRecPoints->GetEntriesFast() )) ; fRecPointsInRun += emcRecPoints->GetEntriesFast() ; fRecPointsInRun += cpvRecPoints->GetEntriesFast() ; @@ -981,14 +1102,14 @@ void AliPHOSClusterizerv1::PrintRecPoints(Option_t * option) printf("\n CPV clusters \n") ; printf("Index Ene(MeV) Module X Y Z \n") ; for (index = 0 ; index < cpvRecPoints->GetEntries() ; index++) { - AliPHOSRecPoint * rp = (AliPHOSRecPoint * )cpvRecPoints->At(index) ; - - TVector3 locpos; - rp->GetLocalPosition(locpos); - - printf("\n%6d %8.2f %2d %4.1f %4.1f %4.1f \n", - rp->GetIndexInList(), rp->GetEnergy(), rp->GetPHOSMod(), - locpos.X(), locpos.Y(), locpos.Z()) ; + AliPHOSCpvRecPoint * rp = (AliPHOSCpvRecPoint * )cpvRecPoints->At(index) ; + + TVector3 locpos; + rp->GetLocalPosition(locpos); + + printf("\n%6d %8.2f %2d %4.1f %4.1f %4.1f \n", + rp->GetIndexInList(), rp->GetEnergy(), rp->GetPHOSMod(), + locpos.X(), locpos.Y(), locpos.Z()) ; } } }