LMU/musrsim
LMU/musrsim
5
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

6.10.2011 - Kamil Sedlak

Browse Source 
1) Added "tubsbox2" as a new volume in musrDetectorConstruction
  2) Changed the way how the double hits in M-counter and P-counter
     are treated in musrSimAna - this was a very big change which 
     still needs to be cross-checked.
This commit is contained in:
sedlak 2011-10-06 09:49:51 +00:00
parent 0c995ad9d8
commit 781fb2022c
5 changed files with 123 additions and 58 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

63
musrSimAna/musrAnalysis.cxx
View File

@ -218,6 +218,12 @@ void musrAnalysis::ReadInInputParameters(char* charV1190FileName) {
bool_debugingRequired=true; bool_debugingRequired=true;
debugEventMap.insert(std::pair<int,int>(ieventToDebug_tmp,iLevelToDebug_tmp)); debugEventMap.insert(std::pair<int,int>(ieventToDebug_tmp,iLevelToDebug_tmp));
} }
else if (strcmp(tmpString0,"USE_UNPERFECT_POSITRONS_IN_DOUBLE_HITS")==0) {
musrCounter::bool_ignoreUnperfectPositrons = false;
}
else if (strcmp(tmpString0,"USE_UNPERFECT_MUONS_IN_DOUBLE_HITS")==0) {
musrCounter::bool_ignoreUnperfectMuons = false;
}
else if (strncmp(tmpString0,"musrTH",strlen("musrTH"))==0) { else if (strncmp(tmpString0,"musrTH",strlen("musrTH"))==0) {
// Definition of the histograms - either musrTH1D or musrTH2D // Definition of the histograms - either musrTH1D or musrTH2D
int beginningOfHistoTitle=0, endOfHistoTitle =0; int beginningOfHistoTitle=0, endOfHistoTitle =0;
@ -1007,12 +1013,13 @@ void musrAnalysis::FillHistograms(Int_t iiiEntry) {
Int_t idetP = 0; Int_t idetP = 0;
Int_t idetP_ID = 0; Int_t idetP_ID = 0;
Double_t idetP_edep = 0.; Double_t idetP_edep = 0.;
Bool_t doubleHitM = false; //cks 4.10.2011 Bool_t doubleHitM = false;
Bool_t doubleHitP = false; //cks 4.10.2011 Bool_t doubleHitP = false;
// std::cout<<" FillHistograms: timeBinOfThePreviouslyProcessedHit = "<<timeBinOfThePreviouslyProcessedHit<<std::endl; // std::cout<<" FillHistograms: timeBinOfThePreviouslyProcessedHit = "<<timeBinOfThePreviouslyProcessedHit<<std::endl;
// mCounterHitExistsForThisEventID = mCounter->GetNextGoodHitInThisEvent(eventID,timeBinOfThePreviouslyProcessedHit,0,'M',timeBin0,kEntry,idetM,idetM_ID,idetM_edep,doubleHitM); // mCounterHitExistsForThisEventID = mCounter->GetNextGoodHitInThisEvent(eventID,timeBinOfThePreviouslyProcessedHit,0,'M',timeBin0,kEntry,idetM,idetM_ID,idetM_edep,doubleHitM);
mCounterHitExistsForThisEventID = mCounter->GetNextGoodMuon(eventID,timeBinOfThePreviouslyProcessedHit,timeBin0,kEntry,idetM,idetM_ID,idetM_edep,doubleHitM); // mCounterHitExistsForThisEventID = mCounter->GetNextGoodMuon(eventID,timeBinOfThePreviouslyProcessedHit,timeBin0,kEntry,idetM,idetM_ID,idetM_edep,doubleHitM);
mCounterHitExistsForThisEventID = MuonCounterHit(eventID,timeBinOfThePreviouslyProcessedHit,timeBin0,kEntry,idetM,idetM_ID,idetM_edep);
timeBinOfThePreviouslyProcessedHit = timeBin0; timeBinOfThePreviouslyProcessedHit = timeBin0;
//___________________________________________________________ //___________________________________________________________
@ -1044,7 +1051,7 @@ void musrAnalysis::FillHistograms(Int_t iiiEntry) {
numberOfGoodMuons++; numberOfGoodMuons++;
Long64_t dataBinMin = timeBin0+dataWindowBinMin; Long64_t dataBinMin = timeBin0+dataWindowBinMin;
Long64_t dataBinMax = timeBin0+dataWindowBinMax; Long64_t dataBinMax = timeBin0+dataWindowBinMax;
pCounterHitExistsForThisEventID = PositronCounterHit(eventID,dataBinMin,dataBinMax,timeBin1,timeBin2,posEntry,idetP,idetP_ID,idetP_edep,doubleHitP); pCounterHitExistsForThisEventID = PositronCounterHit(eventID,dataBinMin,dataBinMax,timeBin1,timeBin2,posEntry,idetP,idetP_ID,idetP_edep);
//cDEL if (pCounterHitExistsForThisEventID) std::cout<<" timeBin1-timeBin2 = "<<timeBin1<<"-"<<timeBin2<<"="<<timeBin1-timeBin2 <<std::endl; //cDEL if (pCounterHitExistsForThisEventID) std::cout<<" timeBin1-timeBin2 = "<<timeBin1<<"-"<<timeBin2<<"="<<timeBin1-timeBin2 <<std::endl;
// //
if (pCounterHitExistsForThisEventID&&(posEntry>0)) { if (pCounterHitExistsForThisEventID&&(posEntry>0)) {
@ -1217,7 +1224,8 @@ void musrAnalysis::FillHistograms(Int_t iiiEntry) {
// Check whether there is a good trigger in the next event // Check whether there is a good trigger in the next event
// std::cout<<" FillHistograms: timeBinOfThePreviouslyProcessedHit = "<<timeBinOfThePreviouslyProcessedHit<<std::endl; // std::cout<<" FillHistograms: timeBinOfThePreviouslyProcessedHit = "<<timeBinOfThePreviouslyProcessedHit<<std::endl;
// mCounterHitExistsForThisEventID = mCounter->GetNextGoodHitInThisEvent(eventID,timeBinOfThePreviouslyProcessedHit,0,'M',timeBin0,kEntry,idetM,idetM_ID,idetM_edep,doubleHitM); // mCounterHitExistsForThisEventID = mCounter->GetNextGoodHitInThisEvent(eventID,timeBinOfThePreviouslyProcessedHit,0,'M',timeBin0,kEntry,idetM,idetM_ID,idetM_edep,doubleHitM);
mCounterHitExistsForThisEventID = mCounter->GetNextGoodMuon(eventID,timeBinOfThePreviouslyProcessedHit,timeBin0,kEntry,idetM,idetM_ID,idetM_edep,doubleHitM); // mCounterHitExistsForThisEventID = mCounter->GetNextGoodMuon(eventID,timeBinOfThePreviouslyProcessedHit,timeBin0,kEntry,idetM,idetM_ID,idetM_edep,doubleHitM);
mCounterHitExistsForThisEventID = MuonCounterHit(eventID,timeBinOfThePreviouslyProcessedHit,timeBin0,kEntry,idetM,idetM_ID,idetM_edep);
timeBinOfThePreviouslyProcessedHit = timeBin0; timeBinOfThePreviouslyProcessedHit = timeBin0;
// if (mCounterHitExistsForThisEventID) std::cout<<" YYYYYYYYYYYYYYYYYYY check this : LOOOPING AGAIN"<<std::endl; // if (mCounterHitExistsForThisEventID) std::cout<<" YYYYYYYYYYYYYYYYYYY check this : LOOOPING AGAIN"<<std::endl;
oncePerEvent=false; oncePerEvent=false;
@ -1355,7 +1363,16 @@ Double_t musrAnalysis::PreprocessEvent(Long64_t iEn) {
} }
//================================================================ //================================================================
Bool_t musrAnalysis::PositronCounterHit(Int_t evID, Long64_t dataBinMin, Long64_t dataBinMax, Long64_t& tBin1, Long64_t& tBin2, Int_t& kEntry, Int_t& idetP, Int_t& idetP_ID, Double_t& idetP_edep, Bool_t& doubleHit) { Bool_t musrAnalysis::MuonCounterHit(Int_t evID, Long64_t timeBinMin, Long64_t& timeBin0, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep) {
Bool_t mCounterHitCanditateExists = mCounter->GetNextGoodMuon(evID,timeBinMin,timeBin0,kEntry,idet,idetID,idetEdep);
if (!mCounterHitCanditateExists) return false;
// Check for other muons within the pileup window:
if ( mCounter->CheckForPileupMuons(timeBin0,kEntry) ) return false; // This muon candidate is killed due to a double hit rejection.
return true;
}
//================================================================
Bool_t musrAnalysis::PositronCounterHit(Int_t evID, Long64_t dataBinMin, Long64_t dataBinMax, Long64_t& tBin1, Long64_t& tBin2, Int_t& kEntry, Int_t& idetP, Int_t& idetP_ID, Double_t& idetP_edep) {
if (bool_debugingRequired) { if (bool_debugingRequired) {
if (debugEventMap[eventID]>2) {std::cout<<"DEBUGEVENT:"<<eventID<<"\"PositronCounterHit\": pCounterMap.size()="<<pCounterMap.size()<<std::endl;} if (debugEventMap[eventID]>2) {std::cout<<"DEBUGEVENT:"<<eventID<<"\"PositronCounterHit\": pCounterMap.size()="<<pCounterMap.size()<<std::endl;}
@ -1364,29 +1381,41 @@ Bool_t musrAnalysis::PositronCounterHit(Int_t evID, Long64_t dataBinMin, Long64_
Bool_t positronHitFound = false; Bool_t positronHitFound = false;
Bool_t goodPositronFound = false; Bool_t goodPositronFound = false;
Int_t positronQuality = 0;
// std::cout<<"Debug 10------------------------------"<<std::endl;
if (musrMode=='D') { if (musrMode=='D') {
// Loop over all positron counters // Loop over all positron counters
for (counterMapType::const_iterator it = pCounterMap.begin(); it!=pCounterMap.end(); ++it) { for (counterMapType::const_iterator it = pCounterMap.begin(); it!=pCounterMap.end(); ++it) {
// std::cout<<" ===POSITRON==="<< pCounterMap.size()<<std::endl; // std::cout<<"Debug 20"<<std::endl;
// Bool_t thereWasHit = (it->second)->GetNextGoodHitInThisEvent(evID,dataBinMin,dataBinMax,'P',tBin1,kEntry,idetP,idetP_ID,idetP_edep,doubleHit); Long64_t dataBinMinimum = dataBinMin;
// Bool_t thereWasHit = (it->second)->GetNextGoodPositron(evID,dataBinMin,dataBinMax,tBin1,tBin2,kEntry,idetP,idetP_ID,idetP_edep,doubleHit); do {
Int_t positronQuality = (it->second)->GetNextGoodPositron(evID,dataBinMin,dataBinMax,tBin1,tBin2,kEntry,idetP,idetP_ID,idetP_edep,doubleHit); // std::cout<<"Debug 30"<<std::endl;
// std::cout<<"000000000000 tBin1="<<tBin1<<" tBin2="<<tBin2<<std::endl; positronQuality = (it->second)->GetNextGoodPositron(evID,dataBinMinimum,dataBinMax,tBin1,tBin2,kEntry,idetP,idetP_ID,idetP_edep);
if (doubleHit) {return false;} // There were two hits in the same positron counter // std::cout<<"Debug 40 positronQuality="<<positronQuality<<std::endl;
dataBinMinimum=tBin1+1; // If a positron is found, one needs to search again through the remaining hits in this counter.
if (positronQuality>0) { if (positronQuality>0) {
if (positronHitFound) { // std::cout<<"Debug 50"<<std::endl;
if (positronHitFound) { // Positron found now but also previously ==> double hit ==> through away this event.
// std::cout<<"Debug 60"<<std::endl;
if (bool_debugingRequired) { if (bool_debugingRequired) {
if (debugEventMap[eventID]>2) {std::cout<<"DEBUGEVENT:"<<eventID<<"\"PositronCounterHit\": Coincidence with other positron candidate in other counter."<<std::endl;} if (debugEventMap[eventID]>2) {std::cout<<"DEBUGEVENT:"<<eventID
<<"\"PositronCounterHit\": Coincidence with other positron candidate in other counter."<<std::endl;}
} }
doubleHit = true;
return false; return false;
} // There were two hits in two different positron counters } // end of "positronHitFound"
// std::cout<<"Debug 70"<<std::endl;
positronHitFound = true; positronHitFound = true;
if (positronQuality>1) goodPositronFound = true; if (positronQuality>1) goodPositronFound = true;
} }
// std::cout<<"Debug 80"<<std::endl;
} while (positronQuality>0);
// std::cout<<"Debug 90"<<std::endl;
} }
// std::cout<<"Debug 100 positronQuality="<<positronQuality<<std::endl;
if (goodPositronFound) return true;
} }
return goodPositronFound; // std::cout<<"Debug 110"<<std::endl;
return false;
} }
//================================================================ //================================================================

3
musrSimAna/musrAnalysis.hh
View File

@ -297,7 +297,8 @@ public :
virtual void PrintHitsInAllCounters(); virtual void PrintHitsInAllCounters();
virtual void InitialiseEvent(); virtual void InitialiseEvent();
virtual Double_t PreprocessEvent(Long64_t iEn); virtual Double_t PreprocessEvent(Long64_t iEn);
virtual Bool_t PositronCounterHit(Int_t evID, Long64_t dataBinMin, Long64_t dataBinMax, Long64_t& tBin1, Long64_t& tBin2, Int_t& kEntry, Int_t& idetP, Int_t& idetP_ID, Double_t& idetP_edep, Bool_t& doubleHit); virtual Bool_t PositronCounterHit(Int_t evID, Long64_t dataBinMin, Long64_t dataBinMax, Long64_t& tBin1, Long64_t& tBin2, Int_t& kEntry, Int_t& idetP, Int_t& idetP_ID, Double_t& idetP_edep);
virtual Bool_t MuonCounterHit(Int_t evID, Long64_t timeBinMin, Long64_t& timeBin0, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep);
void CopySubstring(char* inputChar,int iStart,int iEnd,char* outputChar); void CopySubstring(char* inputChar,int iStart,int iEnd,char* outputChar);
void MyPrintTree(); void MyPrintTree();
void MyPrintConditions(); void MyPrintConditions();

78
musrSimAna/musrCounter.cxx
View File

@ -6,6 +6,9 @@ typedef std::map<int,int> debugEventMapType;
debugEventMapType debugEventMap; debugEventMapType debugEventMap;
Bool_t bool_debugingRequired; Bool_t bool_debugingRequired;
Bool_t musrCounter::bool_ignoreUnperfectMuons = true;
Bool_t musrCounter::bool_ignoreUnperfectPositrons = true;
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
musrCounter::musrCounter(int CHANNEL_NR, char CHANNEL_NAME[200], char CHANNEL_TYPE, float E_THRESH, int TIME_SHIFT) { musrCounter::musrCounter(int CHANNEL_NR, char CHANNEL_NAME[200], char CHANNEL_TYPE, float E_THRESH, int TIME_SHIFT) {
@ -152,12 +155,11 @@ Bool_t musrCounter::IsInCoincidence(Long64_t timeBin, char motherCounter, Bool_t
} }
//================================================================ //================================================================
Bool_t musrCounter::GetNextGoodMuon(Int_t evtID, Long64_t timeBinMin, Long64_t& timeBinOfNextHit, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep, Bool_t& doubleHitFound) { Bool_t musrCounter::GetNextGoodMuon(Int_t evtID, Long64_t timeBinMin, Long64_t& timeBinOfNextHit, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep) {
// This function searches for a good muon, i.e. the muon // This function searches for a good muon, i.e. the muon
// 1) belongs to the currently analysed event // 1) belongs to the currently analysed event
// 2) is in coincidence with all required coincidence detectors // 2) is in coincidence with all required coincidence detectors
// 3) is not in coincidence with veto detectors // 3) is not in coincidence with veto detectors
// 4) is not in coincidence with other muons (more precisely - with hits in m-counter).
// INPUT PARAMETERS: evtID, timeBinMin // INPUT PARAMETERS: evtID, timeBinMin
// OUTPUT PARAMETERS: timeBinOfNextHit // OUTPUT PARAMETERS: timeBinOfNextHit
// //
@ -165,7 +167,6 @@ Bool_t musrCounter::GetNextGoodMuon(Int_t evtID, Long64_t timeBinMin, Long64_t&
// std::cout<<" musrCounter::GetNextGoodMuon timeBinMin="<<timeBinMin<<std::endl; // std::cout<<" musrCounter::GetNextGoodMuon timeBinMin="<<timeBinMin<<std::endl;
if (hitMap.empty()) return false; if (hitMap.empty()) return false;
if (counterType!='M') {std::cout<<"\n!!! FATAL ERROR !!! musrCounter::GetNextGoodMuon: not the muon counter! ==> S T O P !!!\n"; exit(1);} if (counterType!='M') {std::cout<<"\n!!! FATAL ERROR !!! musrCounter::GetNextGoodMuon: not the muon counter! ==> S T O P !!!\n"; exit(1);}
doubleHitFound = false;
for (hitMap_TYPE::iterator it = hitMap.begin(); it != hitMap.end(); ++it) { for (hitMap_TYPE::iterator it = hitMap.begin(); it != hitMap.end(); ++it) {
Long64_t timeBinOfCount_tmp = it->first; Long64_t timeBinOfCount_tmp = it->first;
@ -179,24 +180,20 @@ Bool_t musrCounter::GetNextGoodMuon(Int_t evtID, Long64_t timeBinMin, Long64_t&
// std::cout<<"*** "<<evtID<<" eventNumber="<<eventNumber<<" canditas="<<numberOfMuonCandidates<<std::endl; // std::cout<<"*** "<<evtID<<" eventNumber="<<eventNumber<<" canditas="<<numberOfMuonCandidates<<std::endl;
// Hit candidate was found. Now check its coincidences and vetos // Hit candidate was found. Now check its coincidences and vetos
for (counterMapType::const_iterator itCounter = koincidenceCounterMap.begin(); itCounter!=koincidenceCounterMap.end(); ++itCounter) { for (counterMapType::const_iterator itCounter = koincidenceCounterMap.begin(); itCounter!=koincidenceCounterMap.end(); ++itCounter) {
if (!( (itCounter->second)->IsInCoincidence(timeBinOfCount_tmp,'M') )) goto endOfThisHit; // no coincidence found ==> skip hit if (!( (itCounter->second)->IsInCoincidence(timeBinOfCount_tmp,'M') )) { // no coincidence found ==> skip hit
if (bool_ignoreUnperfectMuons) hitMap.erase(it);
goto endOfThisHit;
}
} }
for (counterMapType::const_iterator itCounter = vetoCounterMap.begin(); itCounter!=vetoCounterMap.end(); ++itCounter) { for (counterMapType::const_iterator itCounter = vetoCounterMap.begin(); itCounter!=vetoCounterMap.end(); ++itCounter) {
if ( (itCounter->second)->IsInCoincidence(timeBinOfCount_tmp,'M') ) goto endOfThisHit; // coincidence with veto found ==> skip hit if ( (itCounter->second)->IsInCoincidence(timeBinOfCount_tmp,'M') ) { // coincidence with veto found ==> skip hit
if (bool_ignoreUnperfectMuons) hitMap.erase(it);
goto endOfThisHit;
}
} }
numberOfMuonCandidatesAfterVK++; numberOfMuonCandidatesAfterVK++;
// Check coincidences with other hits in the M counter
// it is expected that there is only one M counter
if (this->IsInCoincidence(timeBinOfCount_tmp,'M',true) ) { // coincidence with another M-counter hit ==> skip hit
doubleHitFound = true;
// std::cout<<"UJGUR double hit found"<<std::endl;
goto endOfThisHit;
}
numberOfMuonCandidatesAfterVKandDoubleHitRemoval++;
kEntry = (it->second)->eventEntry; kEntry = (it->second)->eventEntry;
idet = (it->second)->det_i; idet = (it->second)->det_i;
idetID = (it->second)->det_id; idetID = (it->second)->det_id;
@ -211,7 +208,39 @@ Bool_t musrCounter::GetNextGoodMuon(Int_t evtID, Long64_t timeBinMin, Long64_t&
} }
//================================================================ //================================================================
Int_t musrCounter::GetNextGoodPositron(Int_t evtID, Long64_t timeBinMin, Long64_t timeBinMax, Long64_t& timeBinOfNextGoodHit, Long64_t& timeBinOfNextGoodHit_phaseShifted, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep, Bool_t& doubleHitFound) { Bool_t musrCounter::CheckForPileupMuons(Long64_t timeBin0, Int_t kEntry_NN) {
// Check for pileup muons. If double hit in M-counter is found, return true.
Long64_t timeBinMinimum = timeBin0;
for (hitMap_TYPE::iterator it = hitMap.begin(); it != hitMap.end(); ++it) {
Long64_t timeBinOfCount_tmp = it->first;
if (timeBinOfCount_tmp < timeBinMinimum+coincidenceTimeWindowMin_M) continue; // This hit happened too long ago
if (timeBinOfCount_tmp > timeBinMinimum+coincidenceTimeWindowMax_M) break; // This hit happened too late
if ((timeBinOfCount_tmp == timeBinMinimum)&&( ((it->second)->eventEntry) == kEntry_NN)) continue;
// This is the M-counter hit for which we check the pileup -> ignore this hit in this double-hit check.
// We have found the hit, which could be the double hit. If it is required, we now have
// to check the coincidences and anticoincidences of this double-hit muon:
if (!bool_ignoreUnperfectMuons) {
return true; // In this case we do not check for coincidences and anticoincidences ==> double hit found.
}
for (counterMapType::const_iterator itCounter = koincidenceCounterMap.begin(); itCounter!=koincidenceCounterMap.end(); ++itCounter) {
if (!( (itCounter->second)->IsInCoincidence(timeBinOfCount_tmp,'M') )) goto endOfThisHit; // no coincidence found ==> skip hit
}
for (counterMapType::const_iterator itCounter = vetoCounterMap.begin(); itCounter!=vetoCounterMap.end(); ++itCounter) {
if ( (itCounter->second)->IsInCoincidence(timeBinOfCount_tmp,'M') ) goto endOfThisHit; // coincidence with veto found ==> skip hit
}
// The double hit was found (the pileup muon fulfils all veto and coincidence requirements) ==> end of the search
return true;
endOfThisHit:
;
}
numberOfMuonCandidatesAfterVKandDoubleHitRemoval++;
return false;
}
//================================================================
Int_t musrCounter::GetNextGoodPositron(Int_t evtID, Long64_t timeBinMin, Long64_t timeBinMax, Long64_t& timeBinOfNextGoodHit, Long64_t& timeBinOfNextGoodHit_phaseShifted, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep) {
// INPUT PARAMETERS: evtID, timeBinMin // INPUT PARAMETERS: evtID, timeBinMin
// OUTPUT PARAMETERS: timeBinOfNextGoodHit // OUTPUT PARAMETERS: timeBinOfNextGoodHit
// positronQuality = 0 ... no positron candidate found // positronQuality = 0 ... no positron candidate found
@ -222,7 +251,6 @@ Int_t musrCounter::GetNextGoodPositron(Int_t evtID, Long64_t timeBinMin, Long64_
if (bool_debugingRequired) {if (debugEventMap[evtID]>4) myPrintThisCounter(evtID);} if (bool_debugingRequired) {if (debugEventMap[evtID]>4) myPrintThisCounter(evtID);}
if (hitMap.empty()) return 0; if (hitMap.empty()) return 0;
if (counterType!='P') {std::cout<<"\n!!! FATAL ERROR !!! musrCounter::GetNextGoodPositron: not the positron counter! ==> S T O P !!!\n"; exit(1);} if (counterType!='P') {std::cout<<"\n!!! FATAL ERROR !!! musrCounter::GetNextGoodPositron: not the positron counter! ==> S T O P !!!\n"; exit(1);}
doubleHitFound = false;
for (hitMap_TYPE::iterator it = hitMap.begin(); it != hitMap.end(); ++it) { for (hitMap_TYPE::iterator it = hitMap.begin(); it != hitMap.end(); ++it) {
Int_t eventNumber = (it->second)->eventIDnumber; Int_t eventNumber = (it->second)->eventIDnumber;
Long64_t timeBinOfCount_tmp = it->first; Long64_t timeBinOfCount_tmp = it->first;
@ -243,6 +271,7 @@ Int_t musrCounter::GetNextGoodPositron(Int_t evtID, Long64_t timeBinMin, Long64_
if (bool_debugingRequired) { if (bool_debugingRequired) {
if (debugEventMap[evtID]>3) {std::cout<<"DEBUGEVENT:"<<evtID<<"\"GetNextGoodPositron\": Coincidence required but not found"<<std::endl;} if (debugEventMap[evtID]>3) {std::cout<<"DEBUGEVENT:"<<evtID<<"\"GetNextGoodPositron\": Coincidence required but not found"<<std::endl;}
} }
if (bool_ignoreUnperfectPositrons) hitMap.erase(it); // no coincidence found ==> remove the candidate.
goto endOfThisHit; // no coincidence found ==> skip hit goto endOfThisHit; // no coincidence found ==> skip hit
} }
} }
@ -251,22 +280,10 @@ Int_t musrCounter::GetNextGoodPositron(Int_t evtID, Long64_t timeBinMin, Long64_
if (bool_debugingRequired) { if (bool_debugingRequired) {
if (debugEventMap[evtID]>3) {std::cout<<"DEBUGEVENT:"<<evtID<<"\"GetNextGoodPositron\": Coincidence vith veto detector found"<<std::endl;} if (debugEventMap[evtID]>3) {std::cout<<"DEBUGEVENT:"<<evtID<<"\"GetNextGoodPositron\": Coincidence vith veto detector found"<<std::endl;}
} }
if (bool_ignoreUnperfectPositrons) hitMap.erase(it); // coincidence with veto found ==> remove the candidate.
goto endOfThisHit; // coincidence with veto found ==> skip hit goto endOfThisHit; // coincidence with veto found ==> skip hit
} }
} }
// Check coincidences with other P counters
// Coincidences with other P counters must be checked in musrAnalysis.cxx
//
// ADD HERE THE CHECK THAT THE POSITRON IS FOUND WITHIN THE DATA INTERVAL, i.e. within <timeBinMin,timeBinMax>.
// AND ALSO CHECK THAT THIS IS TRUE FOR THE COINCIDENCE WITH OTHER POSITRON COUNTERS!!!
if (this->IsInCoincidence(timeBinOfCount_tmp,'P',true,timeBinMin,timeBinMax) ) { // coincidence with another P-counter hit ==> skip hit
if (bool_debugingRequired) {
if (debugEventMap[evtID]>3) {std::cout<<"DEBUGEVENT:"<<evtID<<"\"GetNextGoodPositron\": Coincidence with other positron candidate in the same counter."<<std::endl;}
}
doubleHitFound = true;
// std::cout<<"tttttttttttttttttttt doubleHitN="<<doubleHitN++<<std::endl;
goto endOfThisHit;
}
kEntry = (it->second)->eventEntry; kEntry = (it->second)->eventEntry;
idet = (it->second)->det_i; idet = (it->second)->det_i;
@ -277,7 +294,6 @@ Int_t musrCounter::GetNextGoodPositron(Int_t evtID, Long64_t timeBinMin, Long64_
if (bool_debugingRequired) { if (bool_debugingRequired) {
if (debugEventMap[evtID]>3) {std::cout<<"DEBUGEVENT:"<<evtID<<"\"GetNextGoodPositron\": Good positron candidate found in this counter."<<std::endl;} if (debugEventMap[evtID]>3) {std::cout<<"DEBUGEVENT:"<<evtID<<"\"GetNextGoodPositron\": Good positron candidate found in this counter."<<std::endl;}
} }
//cDEL std::cout<<"timeBinOfNextGoodHit ="<<timeBinOfNextGoodHit<<" timeBinOfNextGoodHit_phaseShifted = "<<timeBinOfNextGoodHit_phaseShifted<<std::endl;
return 2; return 2;
endOfThisHit: endOfThisHit:

9
musrSimAna/musrCounter.hh
View File

@ -59,8 +59,9 @@ class musrCounter {
void RemoveHitsInCounter(Long64_t timeBinLimit); void RemoveHitsInCounter(Long64_t timeBinLimit);
void RewindHitsInCounter(Long64_t timeBinsToRewind); void RewindHitsInCounter(Long64_t timeBinsToRewind);
Bool_t IsInCoincidence(Long64_t timeBin, char motherCounter, Bool_t ignoreHitsAtBinZero=false, Long64_t timeBinMinimum=-123456789, Long64_t timeBinMaximum=-123456789); Bool_t IsInCoincidence(Long64_t timeBin, char motherCounter, Bool_t ignoreHitsAtBinZero=false, Long64_t timeBinMinimum=-123456789, Long64_t timeBinMaximum=-123456789);
Bool_t GetNextGoodMuon(Int_t evtID, Long64_t timeBinMin, Long64_t& timeBinOfNextHit, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep, Bool_t& doubleHitFound); Bool_t GetNextGoodMuon(Int_t evtID, Long64_t timeBinMin, Long64_t& timeBinOfNextHit, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep);
Int_t GetNextGoodPositron(Int_t evtID, Long64_t timeBinMin, Long64_t timeBinMax, Long64_t& timeBinOfNextGoodHit, Long64_t& timeBinOfNextGoodHit_phaseShifted, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep, Bool_t& doubleHitFound); Bool_t CheckForPileupMuons(Long64_t timeBin0, Int_t kEntry_NN);
Int_t GetNextGoodPositron(Int_t evtID, Long64_t timeBinMin, Long64_t timeBinMax, Long64_t& timeBinOfNextGoodHit, Long64_t& timeBinOfNextGoodHit_phaseShifted, Int_t& kEntry, Int_t& idet, Int_t& idetID, Double_t& idetEdep);
void myPrintThisCounter(Int_t evtID, Int_t detail=2); void myPrintThisCounter(Int_t evtID, Int_t detail=2);
Long64_t GetNumberOfMuonCandidates(){return numberOfMuonCandidates;} Long64_t GetNumberOfMuonCandidates(){return numberOfMuonCandidates;}
Long64_t GetNumberOfMuonCandidatesAfterVK(){return numberOfMuonCandidatesAfterVK;} Long64_t GetNumberOfMuonCandidatesAfterVK(){return numberOfMuonCandidatesAfterVK;}
@ -92,6 +93,10 @@ class musrCounter {
Long64_t numberOfMuonCandidates; Long64_t numberOfMuonCandidates;
Long64_t numberOfMuonCandidatesAfterVK; Long64_t numberOfMuonCandidatesAfterVK;
Long64_t numberOfMuonCandidatesAfterVKandDoubleHitRemoval; Long64_t numberOfMuonCandidatesAfterVKandDoubleHitRemoval;
public:
static Bool_t bool_ignoreUnperfectMuons;
static Bool_t bool_ignoreUnperfectPositrons;
}; };
#endif #endif

14
src/musrDetectorConstruction.cc
View File

@ -477,6 +477,19 @@ G4VPhysicalVolume* musrDetectorConstruction::Construct() {
G4Tubs* solidOuterDetTube = new G4Tubs("SolidOuterDetTube",0.*mm,x2*mm,x3*mm,x4*deg,x5*deg); G4Tubs* solidOuterDetTube = new G4Tubs("SolidOuterDetTube",0.*mm,x2*mm,x3*mm,x4*deg,x5*deg);
solid = new G4SubtractionSolid(solidName, solidOuterDetTube, solidInnerDetBox); solid = new G4SubtractionSolid(solidName, solidOuterDetTube, solidInnerDetBox);
} }
else if (strcmp(tmpString2,"tubsbox2")==0){
// Create a tube, from which center a box is cut out. More general box and tube than in "tubsbox".
// x1,x2,x3=box half-widths; x4,x5,x6,x7,x8 define the tube.
sscanf(&line[0],"%*s %*s %*s %s %lf %lf %lf %lf %lf %lf %lf %lf %s",
name,&x1,&x2,&x3,&x4,&x5,&x6,&x7,&x8,material);
sscanf(&line[0],"%*s %*s %*s %*s %*g %*g %*g %*g %*g %*g %*g %*g %*s %lf %lf %lf %s %s %s %d %s",
&posx,&posy,&posz,mothersName,rotMatrix,sensitiveDet,&volumeID,actualFieldName);
solidName+=name;
G4double roundingErr=0.01*mm; // to avoid some displaying problems of the subtracted volumes
G4Box* solidInnerDetBox2 = new G4Box("SolidInnerDetBox2",x1*mm,x2*mm,x3*mm+roundingErr);
G4Tubs* solidOuterDetTube2 = new G4Tubs("SolidOuterDetTube2",x4*mm,x5*mm,x6*mm,x7*deg,x8*deg);
solid = new G4SubtractionSolid(solidName, solidOuterDetTube2, solidInnerDetBox2);
}
else if (strcmp(tmpString2,"GPSforward")==0){ else if (strcmp(tmpString2,"GPSforward")==0){
// Create a box, from which a cone is cut out: // Create a box, from which a cone is cut out:
// x1, x2, x3 = box half-widths; // x1, x2, x3 = box half-widths;
@ -1547,6 +1560,7 @@ void musrDetectorConstruction::DefineMaterials()
new G4Material("ArgonGas", z= 18., a= 39.95*g/mole, density= 0.00000000001*mg/cm3); new G4Material("ArgonGas", z= 18., a= 39.95*g/mole, density= 0.00000000001*mg/cm3);
new G4Material("HeliumGas5mbar",z=2., a=4.002602*g/mole, density= 0.00000088132*g/cm3);
if (musrParameters::boolG4OpticalPhotons) { if (musrParameters::boolG4OpticalPhotons) {