// ESAF : Euso Simulation and Analysis Framework
// $Id: RecoEvent.cc,v 1.18 2005/03/25 18:32:40 naumov Exp $
// Marco Pallavicini created Oct, 16 2003
#include <iostream>
#include "RecoEvent.hh"
#include "RecoModuleData.hh"
#include "RecoCellInfo.hh"
#include "RecoPixelData.hh"
#include "RecoAnalogData.hh"
#include "RecoPhotoElectronData.hh"
#include "TMath.h"
ClassImp(RecoEvent)
// ctor
RecoEvent::RecoEvent(): EsafMsgSource(){
}
// dtor
RecoEvent::~RecoEvent() {
for( size_t i(0); i<fPixelData.size(); i++ ) {
delete fPixelData[i];
}
fPixelData.clear();
for( size_t i(0); i<fPhotoElectronData.size(); i++ ) {
delete fPhotoElectronData[i];
}
fPhotoElectronData.clear();
for( size_t i(0); i<fAnalogData.size(); i++ ) {
delete fAnalogData[i];
}
fAnalogData.clear();
/*map<string, const RecoModuleData*>::iterator it1;
for( it1=fModuleData.begin(); it1!=fModuleData.end(); it1++ ) {
delete (*it1).second;
}*/
fModuleData.clear();
map<Int_t, RecoCellInfo*>::iterator it2;
for( it2=fCellData.begin(); it2!=fCellData.end(); it2++ ) {
delete (*it2).second;
}
fModuleData.clear();
Msg(EsafMsg::Info) << "Event " << GetHeader().GetNum() << " processed and destroyed " << MsgDispatch;
}
// get a specified module data
const RecoModuleData* RecoEvent::GetModuleData( const string& name) {
if ( fModuleData.count(name) > 0 ) {
return fModuleData[name];
} else {
return NULL;
}
}
// put a module data in the map
void RecoEvent::PutRecoModuleData( const string& name, const RecoModuleData* md ) {
fModuleData[name] = md;
}
// return info object of a given macrocell id
RecoCellInfo *RecoEvent::GetRecoCellInfo( Int_t cell_id ) {
if ( fCellData.count( cell_id ) == 0 )
return NULL;
return fCellData[cell_id];
}
// return info of a given photoelectron
RecoPhotoElectronData *RecoEvent::GetRecoPhotoElectronData( Int_t ph_id ) {
if ( ph_id >= 0 && ph_id < fHeader.GetTruePhotoElectrons() ) {
return fPhotoElectronData[ph_id];
} else {
return NULL;
}
}
// return info of a given pixel
RecoPixelData *RecoEvent::GetRecoPixelData( Int_t px_id ) {
if ( px_id >= 0 && px_id < fHeader.GetNumActiveFee() ) {
return fPixelData[px_id];
} else {
return NULL;
}
}
// put RecoCellInfo data
void RecoEvent::PutRecoCellInfo( Int_t cell_id, RecoCellInfo *cell_info ) {
if ( fCellData.count (cell_id) == 0) {
fCellData.insert(map<Int_t, RecoCellInfo*>::
value_type(cell_id, cell_info));
}
else
fCellData[cell_id]=cell_info;
}
// return id of most populous macrocell
Int_t RecoEvent::GetMainMacroCell(Int_t nhits_min=1) {
map<Int_t, Int_t> mcellpx;
for (Int_t i=0; i<(Int_t)fPixelData.size(); i++) {
Int_t cell_id = fPixelData[i]->GetMacroCellId();
if (fPixelData[i]->GetCounts()>=nhits_min) {
if ( mcellpx.count(cell_id) == 0 ) {
mcellpx.insert(map<Int_t, Int_t>::
value_type(cell_id, 1));
} else {
mcellpx[cell_id]++;
}
}
}
map<Int_t, Int_t>::iterator it = mcellpx.begin();
Int_t mmc = (*it).first;
Int_t max = (*it).second;
it++;
for ( ; it != mcellpx.end(); it++) {
if ( (*it).second > max ) {
max = (*it).second;
mmc = (*it).first;
}
}
return mmc;
}
// methods of sorting
void RecoEvent::SortPhGtu() {
Int_t n = fPhotoElectronData.size();
for( Int_t gap=n/2; 0<gap; gap/=2 ) {
for( Int_t i=gap; i<n; i++ ) {
for( Int_t j=i-gap; 0<=j; j-=gap ) {
if( fPhotoElectronData[j+gap]->GetGtu() < fPhotoElectronData[j]->GetGtu() ) {
RecoPhotoElectronData *temp = fPhotoElectronData[j];
fPhotoElectronData[j] = fPhotoElectronData[j+gap];
fPhotoElectronData[j+gap] = temp;
}
}
}
}
}
// Find first GTU in the event
Int_t RecoEvent::GetStartGtu(){
Int_t sGtu = fPixelData[0]->GetGtu();
for(size_t i(1); i<fPixelData.size(); i++) {
if ( fPixelData[i]->GetGtu()<sGtu ) sGtu = fPixelData[i]->GetGtu();
}
return sGtu;
}
// Find last GTU in the event
Int_t RecoEvent::GetEndGtu() {
Int_t eGtu = fPixelData[0]->GetGtu();
for(size_t i(1); i<fPixelData.size(); i++) {
if ( fPixelData[i]->GetGtu()>eGtu ) eGtu = fPixelData[i]->GetGtu();
}
return eGtu;
}
// Give the number of photoelectrons in macrocell cell_id
Int_t RecoEvent::GetNumPh( Int_t cell_id ) {
Int_t count = 0;
for( Int_t i=0; i<(Int_t)fPhotoElectronData.size(); i++ ) {
if ( fPhotoElectronData[i]->GetMacroCell() == cell_id )
count++;
}
return count;
}
// Give the number of active pixels in macrocell cell_id
Int_t RecoEvent::GetNumPx( Int_t cell_id, Int_t hits_minimum = 1 ) {
Int_t count = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hits_minimum) {
count++;
}
}
return count;
}
// find min and max theta of pixel FOV in macrocell cell_id
Float_t RecoEvent::GetThetaMin( Int_t cell_id, Int_t hm = 1 ) {
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
Float_t thmin = fPixelData[dummy[0]]->GetTheta();
for ( Int_t i=1; i<j; i++ ) {
if( fPixelData[dummy[i]]->GetTheta() < thmin )
thmin = fPixelData[dummy[i]]->GetTheta();
}
return thmin;
}
Float_t RecoEvent::GetThetaMax( Int_t cell_id, Int_t hm = 1 ){
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
Float_t thmax = fPixelData[dummy[0]]->GetTheta();
for ( Int_t i=1; i<j; i++ ) {
if( fPixelData[dummy[i]]->GetTheta() > thmax )
thmax = fPixelData[dummy[i]]->GetTheta();
}
return thmax;
}
Float_t RecoEvent::GetThetaRange( Int_t cell_id, Int_t hm = 1 ) {
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
if ( j==0 ) {
hm--;
for( ; j>0 || hm==1; hm-- ) {
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
}
}
if (j==0) return 0.0001;
Float_t thmin = fPixelData[dummy[0]]->GetTheta();
Float_t thmax = fPixelData[dummy[0]]->GetTheta();
Float_t thSigmaMin = fPixelData[dummy[0]]->GetThetaSigma();
Float_t thSigmaMax = fPixelData[dummy[0]]->GetThetaSigma();
for ( Int_t i=1; i<j; i++ ) {
if( fPixelData[dummy[i]]->GetTheta() < thmin ) {
thmin = fPixelData[dummy[i]]->GetTheta();
thSigmaMin = fPixelData[dummy[i]]->GetThetaSigma();
}
if( fPixelData[dummy[i]]->GetTheta() > thmax ) {
thmax = fPixelData[dummy[i]]->GetTheta();
thSigmaMax = fPixelData[dummy[i]]->GetThetaSigma();
}
}
if (TMath::Abs(thmax-thmin)!=0) return TMath::Abs(thmax-thmin);
else return (thSigmaMin+thSigmaMax)/2.;
}
// find min and max phi pos of photoelectrons in macrocell cell_id
Float_t RecoEvent::GetPhiMin( Int_t cell_id, Int_t hm = 1 ) {
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
Float_t phimin = fPixelData[dummy[0]]->GetPhi();
for ( Int_t i=1; i<j; i++ ) {
if( fPixelData[dummy[i]]->GetPhi() < phimin )
phimin = fPixelData[dummy[i]]->GetPhi();
}
return phimin;
}
Float_t RecoEvent::GetPhiMax( Int_t cell_id, Int_t hm = 1 ) {
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
Float_t phimax = fPixelData[dummy[0]]->GetPhi();
for ( Int_t i=1; i<j; i++ ) {
if( fPixelData[dummy[i]]->GetPhi() > phimax )
phimax = fPixelData[dummy[i]]->GetPhi();
}
return phimax;
}
Float_t RecoEvent::GetPhiRange( Int_t cell_id, Int_t hm = 1 ) {
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
if ( j==0 ) {
hm--;
for( ; j>0 || hm==1; hm-- ) {
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
}
}
if (j==0) return 0.0001;
Float_t phimin = fPixelData[dummy[0]]->GetPhi();
Float_t phimax = fPixelData[dummy[0]]->GetPhi();
Float_t phiSigmaMin = fPixelData[dummy[0]]->GetPhiSigma();
Float_t phiSigmaMax = fPixelData[dummy[0]]->GetPhiSigma();
for ( Int_t i=1; i<j; i++ ) {
if( fPixelData[dummy[i]]->GetPhi() < phimin ) {
phimin = fPixelData[dummy[i]]->GetPhi();
phiSigmaMin = fPixelData[dummy[i]]->GetPhiSigma();
}
if( fPixelData[dummy[i]]->GetPhi() > phimax ) {
phimax = fPixelData[dummy[i]]->GetPhi();
phiSigmaMax = fPixelData[dummy[i]]->GetPhiSigma();
}
}
if (TMath::Abs(phimax-phimin)!=0) return TMath::Abs(phimax-phimin);
else return (phiSigmaMin+phiSigmaMax)/2;
}
// find first and last gtu of pixels in macrocell cell_id
Int_t RecoEvent::GetFirstGtu( Int_t cell_id, Int_t hm = 1 ) {
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
Int_t tmin = fPixelData[dummy[0]]->GetGtu();
for ( Int_t i=1; i<j; i++ ) {
if( fPixelData[dummy[i]]->GetGtu() < tmin )
tmin = fPixelData[dummy[i]]->GetGtu();
}
return tmin;
}
Int_t RecoEvent::GetLastGtu( Int_t cell_id, Int_t hm = 1 ) {
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
Int_t tmax = fPixelData[dummy[0]]->GetGtu();
for ( Int_t i=1; i<j; i++ ) {
if( fPixelData[dummy[i]]->GetGtu() > tmax )
tmax = fPixelData[dummy[i]]->GetGtu();
}
return tmax;
}
Int_t RecoEvent::GetGtuRange( Int_t cell_id, Int_t hm = 1 ) {
Int_t dummy[GetNumPx(cell_id)];
Int_t j = 0;
for( Int_t i=0; i<(Int_t)fPixelData.size(); i++ ) {
if ( fPixelData[i]->GetMacroCellId() == cell_id && fPixelData[i]->GetCounts() >= hm ) {
dummy[j] = i;
j++;
}
}
Int_t tmin = fPixelData[dummy[0]]->GetGtu();
Int_t tmax = fPixelData[dummy[0]]->GetGtu();
for ( Int_t i=1; i<j; i++ ) {
if( (fPixelData[dummy[i]]->GetGtu()) < tmin )
tmin = fPixelData[dummy[i]]->GetGtu();
if( (fPixelData[dummy[i]]->GetGtu()) > tmax )
tmax = fPixelData[dummy[i]]->GetGtu();
}
return tmax-tmin+1;
}