Skip to content
Snippets Groups Projects
Commit 3d559ad4 authored by Adrian A. Weber's avatar Adrian A. Weber Committed by Pierre-Alain Loizeau
Browse files

add ICD correction to hitpriducer of RICH and add functionality to Macro.

parent bf6a2a95
No related branches found
No related tags found
1 merge request!826add ICD correction to hitpriducer of RICH and add functionality to Macro.
Pipeline #17528 passed
This diff is collapsed.
......@@ -30,6 +30,7 @@ CbmRichMCbmHitProducer::CbmRichMCbmHitProducer()
,
//fMappingFile("mRICH_Mapping_vert_20190318.geo")
fMappingFile("mRICH_Mapping_vert_20190318_elView.geo")
, fICD_offset_read()
{
}
......@@ -60,6 +61,10 @@ InitStatus CbmRichMCbmHitProducer::Init()
InitMapping();
for (auto& a : fICD_offset_read)
a = 0.;
if (fDoICD) read_ICD(fICD_offset_read, 0);
return kSUCCESS;
}
......@@ -121,7 +126,7 @@ void CbmRichMCbmHitProducer::Exec(Option_t* /*option*/)
void CbmRichMCbmHitProducer::ProcessData(CbmEvent* event)
{
if (event != NULL) {
LOG(info) << "CbmRichMCbmHitProducer CbmEvent mode. CbmEvent # " << event->GetNumber();
LOG(debug) << "CbmRichMCbmHitProducer CbmEvent mode. CbmEvent # " << event->GetNumber();
Int_t nofDigis = event->GetNofData(ECbmDataType::kRichDigi);
//LOG(info) << "nofDigis: " << nofDigis;
......@@ -130,8 +135,8 @@ void CbmRichMCbmHitProducer::ProcessData(CbmEvent* event)
Int_t digiIndex = event->GetIndex(ECbmDataType::kRichDigi, iDigi);
ProcessDigi(event, digiIndex);
}
LOG(info) << "nofDigis: " << nofDigis << "\t\t "
<< "nofHits : " << fNofHits;
LOG(debug) << "nofDigis: " << nofDigis << "\t\t "
<< "nofHits : " << fNofHits;
fNofHits = 0;
}
else {
......@@ -176,13 +181,22 @@ void CbmRichMCbmHitProducer::AddHit(CbmEvent* event, TVector3& posHit, const Cbm
Int_t PmtId)
{
Int_t nofHits = fRichHits->GetEntriesFast();
int32_t DiRICH_Addr = digi->GetAddress();
unsigned int addr = (((DiRICH_Addr >> 24) & 0xF) * 18 * 32) + (((DiRICH_Addr >> 20) & 0xF) * 2 * 32)
+ (((DiRICH_Addr >> 16) & 0xF) * 32) + ((DiRICH_Addr & 0xFFFF) - 0x1);
new ((*fRichHits)[nofHits]) CbmRichHit();
CbmRichHit* hit = (CbmRichHit*) fRichHits->At(nofHits);
hit->SetPosition(posHit);
hit->SetDx(fHitError);
hit->SetDy(fHitError);
hit->SetRefId(index);
hit->SetTime(digi->GetTime());
if (fDoICD) { hit->SetTime(digi->GetTime() - fICD_offset_read.at(addr)); }
else {
hit->SetTime(digi->GetTime());
}
hit->SetToT(digi->GetToT());
hit->SetAddress(digi->GetAddress());
hit->SetPmtId(PmtId);
......@@ -273,4 +287,30 @@ bool CbmRichMCbmHitProducer::RestrictToAerogelAccDec2019(Double_t x, Double_t y)
return true;
}
void CbmRichMCbmHitProducer::read_ICD(std::array<Double_t, 2304>& ICD_offsets, unsigned int iteration)
{
std::string line;
std::ifstream icd_file(Form("icd_offset_it_%u.data", iteration));
unsigned int lineCnt = 0;
if (icd_file.is_open()) {
while (getline(icd_file, line)) {
if (line[0] == '#') continue; // just a comment
std::istringstream iss(line);
unsigned int addr = 0;
Double_t value;
if (!(iss >> addr >> value)) {
LOG(info) << "A Problem accured in line " << lineCnt << "\n";
break;
} // error
lineCnt++;
ICD_offsets.at(addr) += value;
}
icd_file.close();
LOG(info) << "Loaded inter channel delay file icd_offset_it_" << iteration << ".data for RICH.\n";
}
else {
LOG(info) << "Unable to open inter channel delay file icd_offset_it_" << iteration << ".data\n";
}
}
ClassImp(CbmRichMCbmHitProducer)
......@@ -103,6 +103,13 @@ public:
*/
void applyToTCut() { fDoToT = true; }
/**
* Apply correction of the inter channel delay (ICD)
*/
void applyICDCorrection(bool val = true) { fDoICD = val; }
/**
* Apply restriction to Mar2019 mRICH Acceptance (for Simulations)
*/
......@@ -119,6 +126,7 @@ private:
TClonesArray* fRichHits; // RICH hits
TClonesArray* fCbmEvents = nullptr; // CbmEvent for time-based simulations
bool fDoToT = false;
bool fDoICD = false;
bool fRestrictToAcc = false;
bool fRestrictToFullAcc = false;
bool fRestrictToAerogelAccDec2019 = false;
......@@ -135,6 +143,8 @@ private:
std::string fMappingFile;
std::array<Double_t, 2304> fICD_offset_read;
void InitMapping();
bool isInToT(const double ToT);
......@@ -154,6 +164,12 @@ private:
void AddHit(CbmEvent* event, TVector3& posHit, const CbmRichDigi* digi, Int_t index, Int_t PmtId);
/**
* function for loading of a created inter channel delay correction file.
*/
void read_ICD(std::array<Double_t, 2304>& offsets, unsigned int iteration);
/**
* \brief Copy constructor.
*/
......
......@@ -167,10 +167,14 @@ void CbmRichMCbmQaRichOnly::InitHistograms()
fHM->Create2<TH2D>("fhRichRingRadiusChi2", "fhRichRingRadiusChi2; Ring Radius [cm];\\Chi^2 ;;Entries", 70, -0.05,
6.95, 101, 0.0, 10.1);
fHM->Create1<TH1D>("fhHitTimeEvent", "fhHitTimeEvent;time [ns];Entries", 400, -100., 300);
// Digis
fHM->Create2<TH2D>("fhDigisInChnl", "fhDigisInChnl;channel;#Digis;", 2304, -0.5, 2303.5, 50, -0.5, 49.5);
fHM->Create2<TH2D>("fhDigisInDiRICH", "fhDigisInDiRICH;DiRICH;#Digis;", 72, -0.5, 71.5, 300, -0.5, 299.5);
fHM->Create2<TH2D>("fhDigisInChnl", "fhDigisInChnl;channel;#Digis;", 2304, -0.5, 2303.5, 500, -0.5, 499.5);
fHM->Create2<TH2D>("fhDigisInDiRICH", "fhDigisInDiRICH;DiRICH;#Digis;", 72, -0.5, 71.5, 3000, -0.5, 2999.5);
//fHM->Create2<TH2D>("fhDigisTimeTot", "fhDigisTimeTot;LE [ns]; ToT [ns];#Digis;", 200, -50., 150., 300, 15.0, 30.0);
fHM->Create2<TH2D>("fhHitsTimeTot", "fhHitsTimeTot;LE [ns]; ToT [ns];#Digis;", 200, -50., 150., 300, 15.0, 30.0);
}
......@@ -187,7 +191,6 @@ void CbmRichMCbmQaRichOnly::Exec(Option_t* /*option*/)
for (int i = 0; i < fDigiMan->GetNofDigis(ECbmModuleId::kRich); i++) {
const CbmRichDigi* richDigi = fDigiMan->Get<CbmRichDigi>(i);
fHM->H1("fhRichDigisToT")->Fill(richDigi->GetToT());
uint16_t addrDiRICH = (richDigi->GetAddress() >> 16) & 0xFFFF;
uint16_t addrChnl = richDigi->GetAddress() & 0xFFFF;
uint16_t dirichNmbr = ((addrDiRICH >> 8) & 0xF) * 18 + ((addrDiRICH >> 4) & 0xF) * 2 + ((addrDiRICH >> 0) & 0xF);
......@@ -238,6 +241,9 @@ void CbmRichMCbmQaRichOnly::Exec(Option_t* /*option*/)
evRichHitIndx.push_back(iRichHit);
CbmRichHit* richHit = static_cast<CbmRichHit*>(fRichHits->At(iRichHit));
fHM->H1("fhHitTimeEvent")->Fill(richHit->GetTime() - ev->GetStartTime());
fHM->H2("fhHitsTimeTot")->Fill(richHit->GetTime() - ev->GetStartTime(), richHit->GetToT());
uint32_t pmtId = (((richHit->GetAddress()) >> 20) & 0xF) + (((richHit->GetAddress()) >> 24) & 0xF) * 9;
pmtHits[pmtId]++;
......@@ -393,11 +399,25 @@ void CbmRichMCbmQaRichOnly::DrawHist()
DrawH2(fHM->H2("fhDigisInChnl"));
}
// {
// fHM->CreateCanvas("DigisTimeTot", "DigisTimeTot", 600, 600);
// DrawH2(fHM->H2("fhDigisTimeTot"));
// }
{
fHM->CreateCanvas("HitsTimeTot", "HitsTimeTot", 600, 600);
DrawH2(fHM->H2("fhHitsTimeTot"));
}
{
fHM->CreateCanvas("DigisInDiRICH", "DigisInDiRICH", 1200, 600);
DrawH2(fHM->H2("fhDigisInDiRICH"));
}
{
fHM->CreateCanvas("HitTimeEvent", "HitTimeEvent", 1200, 1200);
DrawH1(fHM->H1("fhHitTimeEvent"));
}
{
TCanvas* c = fHM->CreateCanvas("rich_Hits", "rich_Hits", 1200, 1200);
......@@ -575,7 +595,8 @@ void CbmRichMCbmQaRichOnly::Finish()
}
ICD_offset.at(i) += ICD_offset_read.at(i);
}
save_ICD(ICD_offset, 0);
if (fGenerateICD) save_ICD(ICD_offset, 0);
//std::cout<<"Address: "<<InterChannel[0].first << std::endl;
//std::cout<<"Tracks: "<< fTofTracks->GetEntriesFast() <<std::endl;
......
......@@ -104,6 +104,12 @@ public:
*/
void SetTriggerRichHits(Int_t val = 0) { fTriggerRichHits = val; }
/**
* activate the generation of new ICD correction iterations
*/
void SetIcdGeneration(bool val = true) { fGenerateICD = val; }
private:
CbmDigiManager* fDigiMan = nullptr;
......@@ -138,6 +144,8 @@ private:
bool fDoWriteHistToFile = true;
bool fDoDrawCanvas = true;
bool fGenerateICD = false;
std::array<Double_t, 2304> ICD_offset_read;
std::array<Double_t, 2304> ICD_offset;
std::array<uint32_t, 2304> ICD_offset_cnt;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment