, fRecalibrateWithClusterTime(0)
, fMaxEvent(0), fDoTrackMatching(kFALSE)
, fSelectCell(kFALSE), fSelectCellMinE(0), fSelectCellMinFrac(0)
+, fRejectBelowThreshold(kFALSE)
, fRemoveLEDEvents(kTRUE),fRemoveExoticEvents(kFALSE)
-, fImportGeometryFromFile(kFALSE), fImportGeometryFilePath("")
+, fImportGeometryFromFile(kTRUE), fImportGeometryFilePath("")
, fOADBSet(kFALSE), fAccessOADB(kTRUE), fOADBFilePath("")
, fCentralityClass(""), fSelectEMCALEvent(0)
, fEMCALEnergyCut(0.), fEMCALNcellsCut (0)
, fSetCellMCLabelFromCluster(0)
, fRemapMCLabelForAODs(0)
+, fInputFromFilter(0)
{
// Constructor
, fRecalibrateWithClusterTime(0)
, fMaxEvent(0), fDoTrackMatching(kFALSE)
, fSelectCell(kFALSE), fSelectCellMinE(0), fSelectCellMinFrac(0)
+, fRejectBelowThreshold(kFALSE)
, fRemoveLEDEvents(kTRUE), fRemoveExoticEvents(kFALSE)
-, fImportGeometryFromFile(kFALSE), fImportGeometryFilePath("")
+, fImportGeometryFromFile(kTRUE), fImportGeometryFilePath("")
, fOADBSet(kFALSE), fAccessOADB(kTRUE), fOADBFilePath("")
, fCentralityClass(""), fSelectEMCALEvent(0)
, fEMCALEnergyCut(0.), fEMCALNcellsCut (0)
, fSetCellMCLabelFromCluster(0)
, fRemapMCLabelForAODs(0)
+, fInputFromFilter(0)
{
// Constructor
if(fEMCALEnergyCut <= 0) return kTRUE; // accept
- Int_t nCluster = InputEvent() -> GetNumberOfCaloClusters();
- AliVCaloCells * caloCell = InputEvent() -> GetEMCALCells();
- Int_t bc = InputEvent() -> GetBunchCrossNumber();
-
+ Int_t nCluster = fEvent -> GetNumberOfCaloClusters();
+ AliVCaloCells * caloCell = fEvent -> GetEMCALCells();
+ Int_t bc = fEvent -> GetBunchCrossNumber();
+
for(Int_t icalo = 0; icalo < nCluster; icalo++)
{
- AliVCluster *clus = (AliVCluster*) (InputEvent()->GetCaloCluster(icalo));
+ AliVCluster *clus = (AliVCluster*) (fEvent->GetCaloCluster(icalo));
if( ( clus->IsEMCAL() ) && ( clus->GetNCells() > fEMCALNcellsCut ) && ( clus->E() > fEMCALEnergyCut ) &&
fRecoUtils->IsGoodCluster(clus,fGeom,caloCell,bc))
if (fDebug > 0)
printf("AliAnalysisTaskEMCALClusterize::AcceptEventEMCAL() - Accept : E %2.2f > %2.2f, nCells %d > %d \n",
- clus->E(), fEMCALEnergyCut, clus->GetNCells(), fEMCALNcellsCut);
+ clus->E(), fEMCALEnergyCut, clus->GetNCells(), fEMCALNcellsCut);
return kTRUE;
}
return kFALSE;
-}
+}
//_______________________________________________
void AliAnalysisTaskEMCALClusterize::AccessOADB()
if(trecal)
{
- TString passTmp = pass;
- if(pass!="pass1" && pass!="pass2") passTmp = "pass2"; // TEMPORARY FIX FOR LHC11a analysis
- TObjArray *trecalpass=(TObjArray*)trecal->FindObject(passTmp);
+ TString passM = pass;
+ if(pass=="spc_calo") passM = "pass1";
+ TObjArray *trecalpass=(TObjArray*)trecal->FindObject(passM);
if(trecalpass)
{
// Also check if the quality of the event is good if not reject it
fEvent = 0x0;
-
- //Process events if there is a high energy cluster
- if(!AcceptEventEMCAL()) return ;
-
+
AliAODInputHandler* aodIH = dynamic_cast<AliAODInputHandler*>((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
Int_t eventN = Entry();
if(aodIH) eventN = aodIH->GetReadEntry();
AODEvent()->GetEMCALCells()->GetNumberOfCells());
}
}
+ else if(fInputFromFilter)
+ {
+ //printf("Get Input From Filtered AOD\n");
+ fEvent = AODEvent();
+ }
else
{
fEvent = InputEvent();
return ;
}
+ //Process events if there is a high energy cluster
+ if(!AcceptEventEMCAL()) { fEvent = 0x0 ; return ; }
+
//-------------------------------------------------------------------------------------
// Reject events if LED was firing, use only for LHC11a data
// Reject event if triggered by exotic cell and remove exotic cells if not triggered
Int_t label = clus->GetLabel();
Int_t label2 = -1 ;
//printf("Org cluster E %f, Time %e, Index = %d, ID %d, MC label %d\n", clus->E(), clus->GetTOF(),i, clus->GetID(),label );
- //printf("Original list of labels for new cluster : \n");
+ //printf("Original list of labels from old cluster : \n");
//for(Int_t imc = 0; imc < clus->GetNLabels(); imc++) printf("\t Label %d\n",clus->GetLabelAt(imc));
+
if (clus->GetNLabels()>=2) label2 = clus->GetLabelAt(1) ;
UShort_t * index = clus->GetCellsAbsId() ;
for(Int_t icell=0; icell < clus->GetNCells(); icell++ )
fCellTime[index[icell]] = clus->GetTOF();
fCellMatchdEta[index[icell]] = clus->GetTrackDz();
fCellMatchdPhi[index[icell]] = clus->GetTrackDx();
- //printf(" %d,", index[icell] );
}
nClustersOrg++;
}
AliVCaloCells *cells = fEvent->GetEMCALCells();
Int_t bc = InputEvent()->GetBunchCrossNumber();
+
for (Int_t icell = 0; icell < cells->GetNumberOfCells(); icell++)
{
// Get cell values, recalibrate and not include bad channels found in analysis, nor cells with too low energy, nor exotic cell
}
Int_t mcLabel = cells->GetMCLabel(icell);
- //printf("AliAnalysisTaksEMCALClusterize::ClusterizeCells() - cell %d, mc label %d",id,mcLabel);
+ //printf("AliAnalysisTaksEMCALClusterize::ClusterizeCells() - cell %d, mc label %d\n",id,mcLabel);
//if(fCellLabels[id]!=mcLabel)printf("mcLabel %d - %d\n",mcLabel,fCellLabels[id]);
if ( fSetCellMCLabelFromCluster == 1 ) mcLabel = fCellLabels[id]; // Older aliroot MC productions
else if( fSetCellMCLabelFromCluster == 0 && fRemapMCLabelForAODs) RemapMCLabelForAODs(mcLabel);
- //printf(", new label %d\n",mcLabel);
+ else mcLabel = -1; // found later
+
+ //printf("\t new label %d\n",mcLabel);
// Create the digit, put a fake primary deposited energy to trick the clusterizer
// when checking the most likely primary
-
- new((*fDigitsArr)[idigit]) AliEMCALDigit( mcLabel, mcLabel, id, amp, time,AliEMCALDigit::kHG,idigit, 0, 0, 1);
+
+ Float_t efrac = cells->GetEFraction(icell);
+
+ //When checking the MC of digits, give weight to cells with embedded signal
+ if (mcLabel > 0 && efrac < 1.e-6) efrac = 1;
+
+ //printf("******* Cell %d, id %d, e %f, fraction %f, MC label %d, used MC label %d\n",icell,id,amp,cells->GetEFraction(icell),cells->GetMCLabel(icell),mcLabel);
+
+ new((*fDigitsArr)[idigit]) AliEMCALDigit( mcLabel, mcLabel, id, amp, time,AliEMCALDigit::kHG,idigit, 0, 0, amp*efrac);
+ // Last parameter should be MC deposited energy, since it is not available, add just the cell amplitude so that
+ // we give more weight to the MC label of the cell with highest energy in the cluster
idigit++;
}
Double_t cellAmplitude = 0;
Double_t cellTime = 0;
Short_t cellNumber = 0;
- Short_t cellMCLabel = 0;
+ Int_t cellMCLabel = 0;
Double_t cellEFrac = 0;
Int_t nClustersOrg = 0;
//Trigger
header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection
- if (esdevent) header->SetFiredTriggerClasses(esdevent->GetFiredTriggerClasses());
- else if (aodevent) header->SetFiredTriggerClasses(aodevent->GetFiredTriggerClasses());
+
+ header->SetFiredTriggerClasses(fEvent->GetFiredTriggerClasses());
header->SetTriggerMask(fEvent->GetTriggerMask());
header->SetTriggerCluster(fEvent->GetTriggerCluster());
header->SetZDCP2Energy(fEvent->GetZDCP2Energy());
header->SetZDCEMEnergy(fEvent->GetZDCEMEnergy(0),fEvent->GetZDCEMEnergy(1));
- Float_t diamxy[2]={fEvent->GetDiamondX(),fEvent->GetDiamondY()};
+ Float_t diamxy[2]={(Float_t)fEvent->GetDiamondX(),(Float_t)fEvent->GetDiamondY()};
Float_t diamcov[3];
fEvent->GetDiamondCovXY(diamcov);
header->SetDiamond(diamxy,diamcov);
new((*fOutputAODBranch)[i]) AliAODCaloCluster(*newCluster);
- if(DebugLevel() > 1 )
- printf("AliAnalysisTaksEMCALClusterize::UserExec() - New cluster %d of %d, energy %f\n",newCluster->GetID(), kNumberOfCaloClusters, newCluster->E());
+ if(DebugLevel() > 1 )
+ printf("AliAnalysisTaksEMCALClusterize::UserExec() - New cluster %d of %d, energy %f, mc label %d \n",newCluster->GetID(), kNumberOfCaloClusters, newCluster->E(), newCluster->GetLabel());
} // cluster loop
else if (pass.Contains("ass3")) return TString("pass3");
else if (pass.Contains("ass4")) return TString("pass4");
else if (pass.Contains("ass5")) return TString("pass5");
+ else if (pass.Contains("LHC11c") && pass.Contains("spc_calo") ) return TString("spc_calo");
+ else if (pass.Contains("calo") || pass.Contains("high_lumi"))
+ {
+ printf("AliAnalysisTaskEMCALClusterize::GetPass() - Path contains <calo> or <high-lumi>, set as <pass1>\n");
+ return TString("pass1");
+ }
// No condition fullfilled, give a default value
printf("AliAnalysisTaskEMCALClusterize::GetPass() - Pass number string not found \n");
if(fMaxEvent <= 0) fMaxEvent = 1000000000;
if(fSelectCellMinE <= 0) fSelectCellMinE = 0.005;
if(fSelectCellMinFrac <= 0) fSelectCellMinFrac = 0.001;
-
+ fRejectBelowThreshold = kFALSE;
+
//Centrality
if(fCentralityClass == "") fCentralityClass = "V0M";
fCentralityBin[0] = clus->fCentralityBin[0];
fCentralityBin[1] = clus->fCentralityBin[1];
}
-
- // Init geometry, I do not like much to do it like this ...
- if(fImportGeometryFromFile && !gGeoManager)
- {
- if (fImportGeometryFilePath == "") fImportGeometryFilePath = "$ALICE_ROOT/OADB/EMCAL/geometry_2011.root" ; // "$ALICE_ROOT/EVE/alice-data/default_geo.root"
- printf("AliAnalysisTaskEMCALClusterize::Init() - Import %s\n",fImportGeometryFilePath.Data());
- TGeoManager::Import(fImportGeometryFilePath) ;
- }
-
+
}
//_______________________________________________________
fClusterizer->SetPar5 (i, fRecParam->GetPar5(i));
fClusterizer->SetPar6 (i, fRecParam->GetPar6(i));
}//end of loop over parameters
-
+ fClusterizer->SetRejectBelowThreshold(fRejectBelowThreshold);//here we set option of unfolding: split or reject energy
fClusterizer->InitClusterUnfolding();
}// to unfold
if (runnumber < 140000) fGeomName = "EMCAL_FIRSTYEARV1";
else if(runnumber < 171000) fGeomName = "EMCAL_COMPLETEV1";
else fGeomName = "EMCAL_COMPLETE12SMV1";
- printf("AliAnalysisTaskEMCALClusterize::InitGeometry() - Set EMCAL geometry name to <%s> for run %d\n",fGeomName.Data(),runnumber);
+ printf("AliAnalysisTaskEMCALClusterize::InitGeometry() - Set EMCAL geometry name to <%s> for run %d\n",
+ fGeomName.Data(),runnumber);
}
fGeom = AliEMCALGeometry::GetInstance(fGeomName);
+ // Init geometry, I do not like much to do it like this ...
+ if(fImportGeometryFromFile && !gGeoManager)
+ {
+ if(fImportGeometryFilePath=="") // If not specified, set location depending on run number
+ {
+ // "$ALICE_ROOT/EVE/alice-data/default_geo.root"
+ if (runnumber < 140000 &&
+ runnumber >= 100000) fImportGeometryFilePath = "$ALICE_ROOT/OADB/EMCAL/geometry_2010.root";
+ if (runnumber >= 140000 &&
+ runnumber < 171000) fImportGeometryFilePath = "$ALICE_ROOT/OADB/EMCAL/geometry_2011.root";
+ else fImportGeometryFilePath = "$ALICE_ROOT/OADB/EMCAL/geometry_2012.root"; // 2012-2013
+ }
+ printf("AliAnalysisTaskEMCALClusterize::InitGeometry() - Import %s\n",fImportGeometryFilePath.Data());
+ TGeoManager::Import(fImportGeometryFilePath) ;
+ }
+
if(fDebug > 0)
{
printf("AliAnalysisTaskEMCALClusterize::InitGeometry(run=%d)",runnumber);
if(accept && !fRecoUtils->IsExoticCell(absID,cells,bc)) totCellE += ecell;
}
- // TString triggerclasses = "";
- // if(esdevent) triggerclasses = esdevent ->GetFiredTriggerClasses();
- // else triggerclasses = ((AliAODEvent*)event)->GetFiredTriggerClasses();
- // //
+ // TString triggerclasses = event->GetFiredTriggerClasses();
// printf("AliAnalysisTaskEMCALClusterize - reject event %d with cluster - reject event with ncells in SM3 %d and SM4 %d\n",(Int_t)Entry(),ncellsSM3, ncellsSM4);
// if(fFillAODFile) AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kFALSE);;
// return;
if(cells->GetAmplitude(icell) > 0.1 && cells->GetCellNumber(icell)/(24*48)==3) ncellsSM3++;
}
- TString triggerclasses = "";
-
- AliESDEvent *esdevent = dynamic_cast<AliESDEvent*>(fEvent);
- if(esdevent) triggerclasses = esdevent ->GetFiredTriggerClasses();
- else triggerclasses = ((AliAODEvent*)fEvent)->GetFiredTriggerClasses();
+ TString triggerclasses = fEvent->GetFiredTriggerClasses();
Int_t ncellcut = 21;
if(triggerclasses.Contains("EMC")) ncellcut = 35;
absIds[ncellsTrue] = digit->GetId();
ratios[ncellsTrue] = recPoint->GetEnergiesList()[c]/digit->GetAmplitude();
-
+
// In case of unfolding, remove digits with unfolded energy too low
if(fSelectCell)
{
if (recPoint->GetEnergiesList()[c] < fSelectCellMinE || ratios[ncellsTrue] < fSelectCellMinFrac)
{
- if(DebugLevel() > 1)
+ if(DebugLevel() > 1)
{
printf("AliAnalysisTaksEMCALClusterize::RecPoints2Clusters() - Too small energy in cell of cluster: cluster cell %f, digit %f\n",
recPoint->GetEnergiesList()[c],digit->GetAmplitude());
}// Select cells
//Recalculate cluster energy and number of cluster cells in case any of the cells was rejected
- ncellsTrue++;
clusterE +=recPoint->GetEnergiesList()[c];
// In case of embedding, fill ratio with amount of signal,
}//Embedding
+ ncellsTrue++;
+
}// cluster cell loop
if (ncellsTrue < 1)
Int_t parentMult = 0;
Int_t *parentList = recPoint->GetParents(parentMult);
clus->SetLabel(parentList, parentMult);
+// printf("Label list : ");
+// for(Int_t ilabel = 0; ilabel < parentMult; ilabel++ ) printf(" %d ",parentList[ilabel]);
+// printf("\n");
}
} // recPoints loop
void AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs(Int_t & label)
{
// MC label for Cells not remapped after ESD filtering, do it here.
-
+
if(label < 0) return ;
AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fEvent) ;
for(Int_t ind = 0; ind < arr->GetEntriesFast(); ind++ )
{
AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
+ if(!particle) continue ;
+
if(label == particle->Label())
{
label = ind;
label = -1;
//printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - Label not found set to -1 \n");
-
+
}
//________________________________________________
// Get the original clusters that contribute to the new cluster, assign the labels of such clusters
// to the new cluster.
// Only approximatedly valid when output are V1 clusters, handle with care
-
- TArrayI clArray(100) ; //Weird if more than a few clusters are in the origin ...
+
+ TArrayI clArray(300) ; //Weird if more than a few clusters are in the origin ...
clArray.Reset();
Int_t nClu = 0;
Int_t nLabTotOrg = 0;
+ Float_t emax = -1;
+ Int_t idMax = -1;
+
+ AliVEvent * event = fEvent;
+
+ //In case of embedding the MC signal is in the added event
+ AliAODInputHandler* aodIH = dynamic_cast<AliAODInputHandler*>((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
+ if(aodIH && aodIH->GetEventToMerge()) //Embedding
+ event = aodIH->GetEventToMerge(); //Get clusters directly from embedded signal
+
//Find the clusters that originally had the cells
for ( Int_t iLoopCell = 0 ; iLoopCell < clus->GetNCells() ; iLoopCell++ )
Bool_t set = kTRUE;
for(Int_t icl =0; icl < nClu; icl++)
{
+ if(((Int_t)clArray.GetAt(icl))==-1 ) continue;
if( idCluster == ((Int_t)clArray.GetAt(icl)) ) set = kFALSE;
// printf("\t \t icell %d Cluster in array %d, IdCluster %d, in array %d, set %d\n",
// iLoopCell, icl, idCluster,((Int_t)clArray.GetAt(icl)),set);
{
clArray.SetAt(idCluster,nClu++);
//printf("******** idCluster %d \n",idCluster);
- nLabTotOrg+=(fEvent->GetCaloCluster(idCluster))->GetNLabels();
+ nLabTotOrg+=(event->GetCaloCluster(idCluster))->GetNLabels();
+
+ //printf("Cluster in array %d, IdCluster %d\n",nClu-1, idCluster);
+
+ //Search highest E cluster
+ AliVCluster * clOrg = event->GetCaloCluster(idCluster);
+ //printf("\t E %f\n",clOrg->E());
+ if(emax < clOrg->E())
+ {
+ emax = clOrg->E();
+ idMax = idCluster;
+ }
}
}
}// cell loop
-
+
if(nClu==0 || nLabTotOrg == 0)
{
- if(clus->E() > 0.25) printf("AliAnalysisTaskEMCALClusterize::SetClustersMCLabelFromOriginalClusters() - Check: N org clusters %d, n tot labels %d, cluster E %f, n cells %d\n",nClu,nLabTotOrg,clus->E(), clus->GetNCells());
+ //if(clus->E() > 0.25) printf("AliAnalysisTaskEMCALClusterize::SetClustersMCLabelFromOriginalClusters() - Check: N org clusters %d, n tot labels %d, cluster E %f, n cells %d\n",nClu,nLabTotOrg,clus->E(), clus->GetNCells());
//for(Int_t icell = 0; icell < clus->GetNCells(); icell++) printf("\t cell %d",clus->GetCellsAbsId()[icell]);
//printf("\n");
}
+ // Put the first in the list the cluster with highest energy
+ if(idMax != ((Int_t)clArray.GetAt(0))) // Max not at first position
+ {
+ Int_t maxIndex = -1;
+ Int_t firstCluster = ((Int_t)clArray.GetAt(0));
+ for ( Int_t iLoopCluster = 0 ; iLoopCluster < nClu ; iLoopCluster++ )
+ {
+ if(idMax == ((Int_t)clArray.GetAt(iLoopCluster))) maxIndex = iLoopCluster;
+ }
+
+ if(firstCluster >=0 && idMax >=0)
+ {
+ clArray.SetAt(idMax,0);
+ clArray.SetAt(firstCluster,maxIndex);
+ }
+ }
+
// Get the labels list in the original clusters, assign all to the new cluster
TArrayI clMCArray(nLabTotOrg) ;
clMCArray.Reset();
Int_t nLabTot = 0;
for ( Int_t iLoopCluster = 0 ; iLoopCluster < nClu ; iLoopCluster++ )
{
- Int_t idCluster = clArray.GetAt(iLoopCluster);
- //printf("\t cl %d \n",idCluster);
- AliVCluster * clOrg = fEvent->GetCaloCluster(idCluster);
+ Int_t idCluster = (Int_t) clArray.GetAt(iLoopCluster);
+ //printf("New Cluster in Array %d, idCluster %d \n",iLoopCluster,idCluster);
+ AliVCluster * clOrg = event->GetCaloCluster(idCluster);
Int_t nLab = clOrg->GetNLabels();
for ( Int_t iLab = 0 ; iLab < nLab ; iLab++ )
clus->SetLabel(clMCArray.GetArray(), nLabTot);
- //printf("Final list of labels for new cluster : \n");
- //for(Int_t imc = 0; imc < clus->GetNLabels(); imc++) printf("\t Label %d\n",clus->GetLabelAt(imc));
-
+// printf("Final list of labels for new cluster : \n");
+// for(Int_t ice = 0; ice < clus->GetNCells() ; ice++)
+// {
+// printf("\t Cell %d ",clus->GetCellsAbsId()[ice]);
+// Int_t label = InputEvent()->GetEMCALCells()->GetCellMCLabel(clus->GetCellsAbsId()[ice]);
+// printf(" org %d ",label);
+// RemapMCLabelForAODs(label);
+// printf(" new %d \n",label);
+// }
+// for(Int_t imc = 0; imc < clus->GetNLabels(); imc++) printf("\t Label %d\n",clus->GetLabelAt(imc));
}