X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=EMCAL%2FAliEMCALReconstructor.cxx;h=63094c31a38c3006e40da083d4ec3b6def5caba6;hb=e2c14c6e829044a38e9a99af5fb883f0867ce68e;hp=a3a4090604a4055bc377e5472585518185bd6131;hpb=0dc91659aeb7ce5476dfef0bc9cd0a3b8a70a9e0;p=u%2Fmrichter%2FAliRoot.git diff --git a/EMCAL/AliEMCALReconstructor.cxx b/EMCAL/AliEMCALReconstructor.cxx index a3a4090604a..63094c31a38 100644 --- a/EMCAL/AliEMCALReconstructor.cxx +++ b/EMCAL/AliEMCALReconstructor.cxx @@ -53,11 +53,14 @@ #include "AliCDBManager.h" #include "AliEMCALGeometry.h" #include "AliEMCAL.h" -#include "AliEMCALHistoUtilities.h" #include "AliESDVZERO.h" #include "AliCDBManager.h" #include "AliRunLoader.h" #include "AliRun.h" +#include "AliEMCALTriggerData.h" +#include "AliEMCALTriggerElectronics.h" +#include "AliEMCALTriggerDCSConfigDB.h" +#include "AliEMCALTriggerDCSConfig.h" ClassImp(AliEMCALReconstructor) @@ -65,9 +68,10 @@ const AliEMCALRecParam* AliEMCALReconstructor::fgkRecParam = 0; // EMCAL rec. p AliEMCALRawUtils* AliEMCALReconstructor::fgRawUtils = 0; // EMCAL raw utilities class AliEMCALClusterizer* AliEMCALReconstructor::fgClusterizer = 0; // EMCAL clusterizer class TClonesArray* AliEMCALReconstructor::fgDigitsArr = 0; // shoud read just once at event +AliEMCALTriggerElectronics* AliEMCALReconstructor::fgTriggerProcessor = 0x0; //____________________________________________________________________________ AliEMCALReconstructor::AliEMCALReconstructor() - : fDebug(kFALSE), fList(0), fGeom(0),fCalibData(0) + : fDebug(kFALSE), fList(0), fGeom(0),fCalibData(0),fPedestalData(0) { // ctor @@ -94,11 +98,29 @@ AliEMCALReconstructor::AliEMCALReconstructor() if(!fCalibData) AliFatal("Calibration parameters not found in CDB!"); + //Get calibration parameters + if(!fPedestalData) + { + AliCDBEntry *entry = (AliCDBEntry*) + AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals"); + if (entry) fPedestalData = (AliCaloCalibPedestal*) entry->GetObject(); + } + + if(!fPedestalData) + AliFatal("Dead map not found in CDB!"); + + //Init the clusterizer with geometry and calibration pointers, avoid doing it twice. - fgClusterizer = new AliEMCALClusterizerv1(fGeom, fCalibData); + fgClusterizer = new AliEMCALClusterizerv1(fGeom, fCalibData,fPedestalData); if(!fGeom) AliFatal(Form("Could not get geometry!")); + AliEMCALTriggerDCSConfigDB* dcsConfigDB = AliEMCALTriggerDCSConfigDB::Instance(); + + const AliEMCALTriggerDCSConfig* dcsConfig = dcsConfigDB->GetTriggerDCSConfig(); + + if (!dcsConfig) AliFatal("No Trigger DCS Configuration from OCDB!"); + fgTriggerProcessor = new AliEMCALTriggerElectronics( dcsConfig ); } //____________________________________________________________________________ @@ -106,15 +128,19 @@ AliEMCALReconstructor::~AliEMCALReconstructor() { // dtor delete fGeom; + delete fgRawUtils; + delete fgClusterizer; + delete fgTriggerProcessor; + AliCodeTimer::Instance()->Print(); } -//____________________________________________________________________________ -void AliEMCALReconstructor::Init() -{ - // Trigger hists - Oct 24, 2007 - fList = AliEMCALHistoUtilities::GetTriggersListOfHists(kTRUE); -} +// //____________________________________________________________________________ +// void AliEMCALReconstructor::Init() +// { +// // Trigger hists - Oct 24, 2007 +// fList = AliEMCALHistoUtilities::GetTriggersListOfHists(kTRUE); +// } //____________________________________________________________________________ void AliEMCALReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const @@ -125,24 +151,57 @@ void AliEMCALReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) // the global tracking. // Works on the current event. - AliCodeTimerAuto("",) + AliCodeTimerAuto("",0) ReadDigitsArrayFromTree(digitsTree); fgClusterizer->InitParameters(); fgClusterizer->SetOutput(clustersTree); + + AliEMCALTriggerData* trgData = new AliEMCALTriggerData(); + + Int_t bufferSize = 32000; + + if (TBranch* triggerBranch = clustersTree->GetBranch("EMTRG")) + triggerBranch->SetAddress(&trgData); + else + clustersTree->Branch("EMTRG","AliEMCALTriggerData",&trgData,bufferSize); + + TClonesArray *trgDigits = new TClonesArray("AliEMCALRawDigit",1000); + TBranch *branchdig = digitsTree->GetBranch("EMTRG"); + if (!branchdig) + { + AliError("Can't get the branch with the EMCAL trigger digits !"); + return; + } + + branchdig->SetAddress(&trgDigits); + branchdig->GetEntry(0); + + //Skip clusterization of LED events + if (GetRecParam()->GetEventSpecie()!=AliRecoParam::kCalib){ + + Int_t v0M[2] = {0,0}; + fgTriggerProcessor->Digits2Trigger(trgDigits, v0M, trgData); - if(fgDigitsArr && fgDigitsArr->GetEntries()) { + + if(fgDigitsArr && fgDigitsArr->GetEntries()) { - fgClusterizer->SetInput(digitsTree); + fgClusterizer->SetInput(digitsTree); - if(Debug()) - fgClusterizer->Digits2Clusters("deb all") ; - else - fgClusterizer->Digits2Clusters(""); + if(Debug()) + fgClusterizer->Digits2Clusters("deb all") ; + else + fgClusterizer->Digits2Clusters(""); - fgClusterizer->Clear(); + fgClusterizer->Clear(); - } + }//digits array exists and has somethind + }//not a LED event + + clustersTree->Fill(); + trgDigits->Delete(); + delete trgDigits; trgDigits = 0x0; + delete trgData; trgData = 0x0; } @@ -157,23 +216,41 @@ void AliEMCALReconstructor::ConvertDigits(AliRawReader* rawReader, TTree* digits rawReader->Reset() ; TClonesArray *digitsArr = new TClonesArray("AliEMCALDigit",200); + TClonesArray *digitsTrg = new TClonesArray("AliEMCALRawDigit", 200); + Int_t bufsize = 32000; digitsTree->Branch("EMCAL", &digitsArr, bufsize); - - //must be done here because, in constructor, option is not yet known - fgRawUtils->SetOption(GetOption()); - - fgRawUtils->SetRawFormatHighLowGainFactor(GetRecParam()->GetHighLowGainFactor()); - fgRawUtils->SetRawFormatOrder(GetRecParam()->GetOrderParameter()); - fgRawUtils->SetRawFormatTau(GetRecParam()->GetTau()); - fgRawUtils->SetNoiseThreshold(GetRecParam()->GetNoiseThreshold()); - fgRawUtils->SetNPedSamples(GetRecParam()->GetNPedSamples()); - - fgRawUtils->Raw2Digits(rawReader,digitsArr); - + digitsTree->Branch("EMTRG", &digitsTrg, bufsize); + + //Skip calibration events do the rest + Bool_t doFit = kTRUE; + if ( !(GetRecParam()->FitLEDEvents()) && GetRecParam()->GetEventSpecie()==AliRecoParam::kCalib) doFit = kFALSE; + if (doFit){ + //must be done here because, in constructor, option is not yet known + fgRawUtils->SetOption(GetOption()); + + fgRawUtils->SetRawFormatHighLowGainFactor(GetRecParam()->GetHighLowGainFactor()); + fgRawUtils->SetRawFormatOrder(GetRecParam()->GetOrderParameter()); + fgRawUtils->SetRawFormatTau(GetRecParam()->GetTau()); + fgRawUtils->SetNoiseThreshold(GetRecParam()->GetNoiseThreshold()); + fgRawUtils->SetNPedSamples(GetRecParam()->GetNPedSamples()); + fgRawUtils->SetRemoveBadChannels(GetRecParam()->GetRemoveBadChannels()); + fgRawUtils->SetFittingAlgorithm(GetRecParam()->GetFittingAlgorithm()); + fgRawUtils->SetFALTROUsage(GetRecParam()->UseFALTRO()); + fgRawUtils->SetTimeMin(GetRecParam()->GetTimeMin()); + fgRawUtils->SetTimeMax(GetRecParam()->GetTimeMax()); + + fgRawUtils->Raw2Digits(rawReader,digitsArr,fPedestalData,digitsTrg); + }//skip calibration event + else{ + AliDebug(1," Calibration Event, skip!"); + } + digitsTree->Fill(); digitsArr->Delete(); + digitsTrg->Delete(); delete digitsArr; + delete digitsTrg; } @@ -188,89 +265,6 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, // printf(" ## AliEMCALReconstructor::FillESD() is started ### \n "); //return; - //###################################################### - //#########Calculate trigger and set trigger info########### - //###################################################### - - AliEMCALTrigger tr; - // tr.SetPatchSize(1); // create 4x4 patches - tr.SetSimulation(kFALSE); // Reconstruction mode - tr.SetDigitsList(fgDigitsArr); - // Get VZERO total multiplicity for jet trigger simulation - // The simulation of jey trigger will be incorrect if no VZERO data - // at ESD - AliESDVZERO* vZero = esd->GetVZEROData(); - if(vZero) { - tr.SetVZER0Multiplicity(vZero->GetMTotV0A() + vZero->GetMTotV0C()); - } - // - tr.Trigger(); - - Float_t maxAmp2x2 = tr.Get2x2MaxAmplitude(); - Float_t maxAmpnxn = tr.GetnxnMaxAmplitude(); - Float_t ampOutOfPatch2x2 = tr.Get2x2AmpOutOfPatch() ; - Float_t ampOutOfPatchnxn = tr.GetnxnAmpOutOfPatch() ; - - Int_t iSM2x2 = tr.Get2x2SuperModule(); - Int_t iSMnxn = tr.GetnxnSuperModule(); - Int_t iModulePhi2x2 = tr.Get2x2ModulePhi(); - Int_t iModulePhinxn = tr.GetnxnModulePhi(); - Int_t iModuleEta2x2 = tr.Get2x2ModuleEta(); - Int_t iModuleEtanxn = tr.GetnxnModuleEta(); - - AliDebug(2, Form("Trigger 2x2 max amp %f, out amp %f, SM %d, iphi %d ieta %d", maxAmp2x2, ampOutOfPatch2x2, iSM2x2,iModulePhi2x2, iModuleEta2x2)); - AliDebug(2, Form("Trigger 4x4 max amp %f , out amp %f, SM %d, iphi %d, ieta %d", maxAmpnxn, ampOutOfPatchnxn, iSMnxn,iModulePhinxn, iModuleEtanxn)); - - TVector3 pos2x2(-1,-1,-1); - TVector3 posnxn(-1,-1,-1); - - Int_t iAbsId2x2 = fGeom->GetAbsCellIdFromCellIndexes( iSM2x2, iModulePhi2x2, iModuleEta2x2) ; // should be changed to Module - Int_t iAbsIdnxn = fGeom->GetAbsCellIdFromCellIndexes( iSMnxn, iModulePhinxn, iModuleEtanxn) ; - fGeom->GetGlobal(iAbsId2x2, pos2x2); - fGeom->GetGlobal(iAbsIdnxn, posnxn); - //printf(" iAbsId2x2 %i iAbsIdnxn %i \n", iAbsId2x2, iAbsIdnxn); - - 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) ; - //printf(" triggerPosition "); - //for(int i=0; i<6; i++) printf(" %i %f : ", i, triggerPosition[i]); - - TArrayF triggerAmplitudes(4); - triggerAmplitudes[0] = maxAmp2x2 ; - triggerAmplitudes[1] = ampOutOfPatch2x2 ; - triggerAmplitudes[2] = maxAmpnxn ; - triggerAmplitudes[3] = ampOutOfPatchnxn ; - //printf("\n triggerAmplitudes "); - //for(int i=0; i<4; i++) printf(" %i %f : ", i, triggerAmplitudes[i]); - //printf("\n"); - //tr.Print(""); - // - // Trigger jet staff - // - if(tr.GetNJetThreshold()>0) { - // Jet phi/eta - Int_t n0 = triggerPosition.GetSize(); - const TH2F *hpatch = tr.GetJetMatrixE(); - triggerPosition.Set(n0 + 2); - for(Int_t i=0; i<2; i++) triggerPosition[n0+i] = hpatch->GetMean(i+1); - // Add jet ampitudes - n0 = triggerAmplitudes.GetSize(); - triggerAmplitudes.Set(n0 + tr.GetNJetThreshold()); - Double_t *ampJet = tr.GetL1JetThresholds(); - for(Int_t i=0; iAddEMCALTriggerPosition(triggerPosition); - esd->AddEMCALTriggerAmplitudes(triggerAmplitudes); - // Fill trigger hists - AliEMCALHistoUtilities::FillTriggersListOfHists(fList,&triggerPosition,&triggerAmplitudes); - //######################################## //##############Fill CaloCells############### //######################################## @@ -292,10 +286,12 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, Float_t energy = 0; for (Int_t idig = 0 ; idig < nDigits ; idig++) { const AliEMCALDigit * dig = (const AliEMCALDigit*)digits->At(idig); - if(dig->GetAmp() > 0 ){ - energy = (static_cast (fgClusterizer))->Calibrate(dig->GetAmp(),dig->GetId()); - emcCells.SetCell(idignew,dig->GetId(),energy, dig->GetTime()); - idignew++; + if(dig->GetAmplitude() > 0 ){ + energy = (static_cast (fgClusterizer))->Calibrate(dig->GetAmplitude(),dig->GetTime(),dig->GetId()); //TimeR or Time? + if(energy > 0){ //Digits tagged as bad (dead, hot, not alive) are set to 0 in calibrate, remove them + emcCells.SetCell(idignew,dig->GetId(),energy, dig->GetTime()); + idignew++; + } } } emcCells.SetNumberOfCells(idignew); @@ -304,15 +300,17 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, //------------------------------------------------------------ //-----------------CLUSTERS----------------------------- //------------------------------------------------------------ + clustersTree->SetBranchStatus("*",0); //disable all branches + clustersTree->SetBranchStatus("EMCALECARP",1); //Enable only the branch we need + TObjArray *clusters = new TObjArray(100); TBranch *branch = clustersTree->GetBranch("EMCALECARP"); branch->SetAddress(&clusters); - clustersTree->GetEvent(0); + branch->GetEntry(0); + //clustersTree->GetEvent(0); Int_t nClusters = clusters->GetEntries(), nClustersNew=0; AliDebug(1,Form("%d clusters",nClusters)); - esd->SetFirstEMCALCluster(esd->GetNumberOfCaloClusters()); // Put after Phos clusters - //###################################################### //#######################TRACK MATCHING############### @@ -336,7 +334,6 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, //######################################## //##############Fill CaloClusters############# //######################################## - esd->SetNumberOfEMCALClusters(nClusters); for (Int_t iClust = 0 ; iClust < nClusters ; iClust++) { const AliEMCALRecPoint * clust = (const AliEMCALRecPoint*)clusters->At(iClust); //if(clust->GetClusterType()== AliESDCaloCluster::kEMCALClusterv1) nRP++; else nPC++; @@ -383,6 +380,10 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, ec->SetClusterType(AliESDCaloCluster::kEMCALClusterv1); ec->SetPosition(xyz); ec->SetE(clust->GetEnergy()); + + //Distance to the nearest bad crystal + ec->SetDistanceToBadChannel(clust->GetDistanceToBadTower()); + ec->SetNCells(newCellMult); //Change type of list from short to ushort UShort_t *newAbsIdList = new UShort_t[newCellMult]; @@ -416,24 +417,18 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, delete [] matchedTrack; - esd->SetNumberOfEMCALClusters(nClustersNew); - //if(nClustersNew != nClusters) - //printf(" ##### nClusters %i -> new %i ##### \n", nClusters, nClustersNew ); - //Fill ESDCaloCluster with PID weights - AliEMCALPID *pid = new AliEMCALPID; - //pid->SetPrintInfo(kTRUE); - pid->SetReconstructor(kTRUE); - pid->RunPID(esd); - delete pid; + AliEMCALPID *pid = new AliEMCALPID; + //pid->SetPrintInfo(kTRUE); + pid->SetReconstructor(kTRUE); + pid->RunPID(esd); + delete pid; - delete digits; - delete clusters; + delete digits; + delete clusters; - // printf(" ## AliEMCALReconstructor::FillESD() is ended : ncl %i -> %i ### \n ",nClusters, nClustersNew); - - //Store EMCAL misalignment matrixes - FillMisalMatrixes(esd) ; + //Store EMCAL misalignment matrixes + FillMisalMatrixes(esd) ; } @@ -442,7 +437,7 @@ void AliEMCALReconstructor::FillMisalMatrixes(AliESDEvent* esd)const{ //Store EMCAL matrixes in ESD Header //Check, if matrixes was already stored - for(Int_t sm = 0 ; sm < 12; sm++){ + for(Int_t sm = 0 ; sm < fGeom->GetNumberOfSuperModules(); sm++){ if(esd->GetEMCALMatrix(sm)!=0) return ; } @@ -454,20 +449,24 @@ void AliEMCALReconstructor::FillMisalMatrixes(AliESDEvent* esd)const{ } //Note, that owner of copied marixes will be header char path[255] ; - TGeoHMatrix * m ; - for(Int_t sm = 0; sm < 12; sm++){ + TGeoHMatrix * m = 0x0; + for(Int_t sm = 0; sm < fGeom->GetNumberOfSuperModules(); sm++){ sprintf(path,"/ALIC_1/XEN1_1/SMOD_%d",sm+1) ; //In Geometry modules numbered 1,2,.,5 if(sm >= 10) sprintf(path,"/ALIC_1/XEN1_1/SM10_%d",sm-10+1) ; - if (gGeoManager->cd(path)){ + if (gGeoManager->CheckPath(path)){ + gGeoManager->cd(path); m = gGeoManager->GetCurrentMatrix() ; +// printf("================================================= \n"); +// printf("AliEMCALReconstructor::FixMisalMatrixes(), sm %d, \n",sm); +// m->Print(""); esd->SetEMCALMatrix(new TGeoHMatrix(*m),sm) ; +// printf("================================================= \n"); } else{ esd->SetEMCALMatrix(NULL,sm) ; } } - } @@ -492,3 +491,4 @@ void AliEMCALReconstructor::ReadDigitsArrayFromTree(TTree *digitsTree) const branch->GetEntry(0); } +