// // ******************************************************************** // * License and Disclaimer * // * * // * The GAMOS software is copyright of the Copyright Holders of * // * the GAMOS Collaboration. It is provided under the terms and * // * conditions of the GAMOS Software License, included in the file * // * LICENSE and available at http://fismed.ciemat.es/GAMOS/license .* // * These include a list of copyright holders. * // * * // * Neither the authors of this software system, nor their employing * // * institutes,nor the agencies providing financial support for this * // * work make any representation or warranty, express or implied, * // * regarding this software system or assume any liability for its * // * use. Please see the license in the file LICENSE and URL above * // * for the full disclaimer and the limitation of liability. * // * * // * This code implementation is the result of the scientific and * // * technical work of the GAMOS collaboration. * // * By using, copying, modifying or distributing the software (or * // * any work based on the software) you agree to acknowledge its * // * use in resulting scientific publications, and indicate your * // * acceptance of all terms of the GAMOS Software license. * // ******************************************************************** // #include #include "SHGetHadronicChannelXSUA.hh" #include "SHUtils.hh" #include "GamosCore/GamosBase/Base/include/GmParameterMgr.hh" #include "GamosCore/GamosBase/Base/include/GmAnalysisMgr.hh" #include "GamosCore/GamosGenerator/include/GmGenerator.hh" #include "GamosCore/GamosUtils/include/GmGenUtils.hh" #include "GamosCore/GamosUtils/include/GmG4Utils.hh" #include "GamosCore/GamosGeometry/include/GmGeometryUtils.hh" #include "GamosCore/GamosUtils/include/GmNumberOfEvent.hh" #include "GamosCore/GamosGenerator/include/GmGeneratorChangeEnergyAndMaterial.hh" #include "G4Material.hh" #include "G4Event.hh" #include "G4EventManager.hh" #include "G4RunManager.hh" #include "G4Box.hh" #include "G4LogicalVolume.hh" #include "G4PVPlacement.hh" #include "G4GeometryManager.hh" #include "G4VUserPrimaryGeneratorAction.hh" #include "G4StateManager.hh" #include "G4ParticleTable.hh" #include "G4Material.hh" #include "G4ElementVector.hh" #include "G4HadronicProcess.hh" #include "G4UserRunAction.hh" #include "G4IonTable.hh" #include //------------------------------------------------------------------ SHGetHadronicChannelXSUA::SHGetHadronicChannelXSUA() { } //------------------------------------------------------------------ void SHGetHadronicChannelXSUA::BeginOfRunAction(const G4Run*) { theCrossSections.clear(); bFirstEventInRun = TRUE; bInitialised = FALSE; thePrimaryParticle = 0; theNEventsInRun = 0; if( !bInitialised ) { bOnlyIons = G4bool(GmParameterMgr::GetInstance()->GetNumericValue(theName+":OnlyIons",0)); bNotGammas = G4bool(GmParameterMgr::GetInstance()->GetNumericValue(theName+":NotGammas",0)); bParticleYields = G4bool(GmParameterMgr::GetInstance()->GetNumericValue(theName+":ParticleYields",0)); bHistos = G4bool(GmParameterMgr::GetInstance()->GetNumericValue(theName+":Histos",1)); if( bHistos ) { const G4VUserPrimaryGeneratorAction* generator = G4RunManager::GetRunManager()->GetUserPrimaryGeneratorAction(); GmGeneratorChangeEnergyAndMaterial* generatorCEM = (GmGeneratorChangeEnergyAndMaterial*)(generator); //--- for use in histogram theMinE = generatorCEM->GetMinE(); theMaxE = generatorCEM->GetMaxE(); theNStepsE = generatorCEM->GetNStepsE(); if( !bLogE ) { theMinE -= 0.5 * generatorCEM->GetStepE(); theMaxE += 0.5 * generatorCEM->GetStepE(); } else { theMinE = log10(generatorCEM->GetMinE()) - 0.5 * generatorCEM->GetStepE(); theMaxE = log10(generatorCEM->GetMaxE()) + 0.5 * generatorCEM->GetStepE(); } theMateIdxMax = G4int(GmParameterMgr::GetInstance()->GetNumericValue(theName+":MaterialIndexMax",10)); theChanIdxMax = G4int(GmParameterMgr::GetInstance()->GetNumericValue(theName+":ChannelIndexMax",200)); theProdIdxMax = G4int(GmParameterMgr::GetInstance()->GetNumericValue(theName+":ProductIndexMax",100)); theTargIdxMax = G4int(GmParameterMgr::GetInstance()->GetNumericValue(theName+":TargetNucleusIndexMax",50)); SetHistoNameAndNumber(theName,theFilters,theClassifier); } } } //------------------------------------------------------------------ void SHGetHadronicChannelXSUA::BeginOfEventAction(const G4Event* evt) { if( !bInitialised ){ G4PrimaryParticle* primParticle = evt->GetPrimaryVertex(0)->GetPrimary(0); SHUtils::StoreHadronicProcesses(primParticle->GetParticleDefinition()->GetParticleName()); if( bHistos ) { const G4VUserPrimaryGeneratorAction* generator = G4RunManager::GetRunManager()->GetUserPrimaryGeneratorAction(); GmGeneratorChangeEnergyAndMaterial* generatorCEM = (GmGeneratorChangeEnergyAndMaterial*)(generator); // create histograms for TOTAL_XS std::vector mates = generatorCEM->GetMaterials(); G4cout << "NMATES " << mates.size() << " NPROC " << SHUtils::GetNHadronicProcesses() << G4endl; for( unsigned int imat = 0; imat < mates.size(); imat++ ) { theMaterials[mates[imat]] = imat; } for( unsigned int imat = 0; imat < mates.size(); imat++ ) { for( G4int iproc = 0; iproc < SHUtils::GetNHadronicProcesses(); iproc++ ){ G4HadronicProcess* hadProc = SHUtils::GetHadronicProcess(iproc); G4String hname = mates[imat]->GetName() + "_" + hadProc->GetProcessName(); // G4cout << " PROCNAME " << hadProc->GetProcessName() << G4endl; theAnaMgr->CreateHisto1D((hname+"_XSmicro").c_str(), theNStepsE, theMinE, theMaxE, theHistoNumber+(imat)*100000+iproc*10000+1); theAnaMgr->CreateHisto1D((hname+"_XSmacro").c_str(), theNStepsE, theMinE, theMaxE, theHistoNumber+(imat)*100000+iproc*10000+2); } } } bInitialised = TRUE; } const G4ParticleDefinition* primary =evt->GetPrimaryVertex(0)->GetPrimary(0)->GetParticleDefinition(); if( thePrimaryParticle != 0 && primary != thePrimaryParticle ) { G4Exception("SHGetHadronicChannelXSUA::BeginOfEventAction", "Warning", JustWarning, G4String("First primary of current event different than previous event : " + primary->GetParticleName() + " <> " + thePrimaryParticle->GetParticleName()).c_str()); } thePrimaryParticle = primary; thePrimaryParticleName = thePrimaryParticle->GetParticleName(); if( bFirstEventInRun ) { InitialiseProcessSecoParticles(); G4PrimaryParticle* particle = evt->GetPrimaryVertex(0)->GetPrimary(0); thePrimaryKineticEnergy = GmGenUtils::GetKineticEnergy(particle->GetMass(), particle->GetMomentum().mag() ); theCurrentMaterial = GmGeometryUtils::GetInstance()->GetTopLV()->GetMaterial(); theCurrentElementVector = const_cast(theCurrentMaterial->GetElementVector()); theCurrentAtomNumDensityVector = const_cast(theCurrentMaterial->GetVecNbOfAtomsPerVolume()); } theNEventsInRun++; } //------------------------------------------------------------------ void SHGetHadronicChannelXSUA::InitialiseProcessSecoParticles() { theProcessSecoParticles.clear(); //-- Get process manager of particle G4ParticleTable* partTable = G4ParticleTable::GetParticleTable(); G4ParticleTable::G4PTblDicIterator* itepar = partTable->GetIterator(); G4ProcessManager* pmanager = 0; itepar->reset(); while( (*itepar)() ){ G4ParticleDefinition* particle = itepar->value(); // G4cout << ":SHGetHadronicChannelXSUA PARTICLE " << particle->GetParticleName() << " lepton " << particle->GetLeptonNumber() << G4endl; if( particle->GetParticleName() == thePrimaryParticleName ){ // G4cout << "SHGetHadronicChannelXSUA: PARTICLE is " << particle->GetParticleName() << G4endl; pmanager = particle->GetProcessManager(); G4ProcessVector* pvect = pmanager->GetProcessList(); int jj, sizproc = pvect->size(); for( jj = 0; jj < sizproc; jj++ ) { G4VProcess* proc = (*pvect)[jj]; if( proc->GetProcessType() != fHadronic ) continue; G4HadronicProcess* prochad = dynamic_cast(proc); if( prochad ) { // G4cout << "SHGetHadronicChannelXSUA::InitialiseProcessSecoParticles proc= " << prochad->GetProcessName() << G4endl; std::map< std::map, G4int > pseco; theProcessSecoParticles[prochad] = pseco; } } } } } //------------------------------------------------------------------ void SHGetHadronicChannelXSUA::UserSteppingAction(const G4Step* aStep) { const G4VProcess* proc = aStep->GetPostStepPoint()->GetProcessDefinedStep(); //----- Only for selected particle const G4Track* aTrack = aStep->GetTrack(); G4Track* aTracknc = const_cast(aTrack); aTracknc->SetTrackStatus( fStopAndKill); // G4cout << " SHGetHadronicChannelXSUA::UserSteppingAction " << aTrack->GetDefinition()->GetParticleName() << " !=? " << theParticleName << G4endl; if( proc->GetProcessType() != fHadronic ) return; if( aTrack->GetDefinition()->GetParticleName() != thePrimaryParticleName ) return; //----- Loop to processes if( bFirstEventInRun ) { // only once per world bFirstEventInRun = FALSE; for( G4int iproc = 0; iproc < SHUtils::GetNHadronicProcesses(); iproc++ ){ G4HadronicProcess* hadProc = SHUtils::GetHadronicProcess(iproc); //--- Get Mean free path aTracknc->SetKineticEnergy( aTracknc->GetVertexKineticEnergy() ); // G4ForceCondition* cond; // G4double mfp = hadProc->GetMeanFreePath(*aTracknc,0.,cond ); //--- Get microscopic cross section const G4DynamicParticle* aParticle = aTracknc->GetDynamicParticle(); G4double cross = 0.; size_t nelm = theCurrentMaterial->GetNumberOfElements(); // G4double aTemp = theCurrentMaterial->GetTemperature(); G4double nAtomsPermm3 = 0; for (size_t iele=0; ieleGetInelasticCrossSectionPerAtom(aParticle,ener,elm); hadProc->GetMicroscopicCrossSection(aParticle,elm,theCurrentMaterial); // G4cout << " theCurrentAtomNumDensityVector " <GetMicroscopicCrossSection(aParticle,elm,aTemp) << G4endl; nAtomsPermm3 += theCurrentAtomNumDensityVector[iele]; } theCrossSections[hadProc] = cross/nAtomsPermm3/barn; G4cout << "@@SHGetHadronicChannelXSUA " << theCurrentMaterial->GetName() << " " << thePrimaryKineticEnergy << " MeV " << hadProc->GetProcessName() << " TOTAL_XS(macro)= " << cross << " mm-1" << " TOTAL_XS(micro)= " << cross/nAtomsPermm3/barn << " barn" << G4endl; if( bHistos ) { G4int materialID = theMaterials[theCurrentMaterial]; GmHisto1* hisXSmicro = theAnaMgr->GetHisto1( theHistoNumber+materialID*100000+iproc*10000+1); hisXSmicro->Fill(thePrimaryKineticEnergy,cross); GmHisto1* hisXSmacro = theAnaMgr->GetHisto1( theHistoNumber+materialID*100000+iproc*10000+2); hisXSmacro->Fill(thePrimaryKineticEnergy,cross/nAtomsPermm3/barn); } } } //---- Get secondary particles const G4TrackVector* secoList = fpSteppingManager->GetSecondary(); //only one step is allowed G4TrackVector::const_iterator ite; std::map secoParticles; // Add target nucleus as first particle in list of secondary particles const G4HadronicProcess* hadProc = (const G4HadronicProcess*)(proc); const G4Nucleus* aNucleus = hadProc->GetTargetNucleus(); // G4String nucName = "[" + GmG4Utils::GetElementNameFromZ(aNucleus->GetZ_asInt() ) + GmG4Utils::itoa(aNucleus->GetA_asInt() ) + " "; // add [ for identification as ion G4String targetNucName = GmG4Utils::GetElementNameFromZ(aNucleus->GetZ_asInt() ) + GmGenUtils::itoa(aNucleus->GetA_asInt() ); theNEvents[targetNucName] ++; targetNucName = "[" + targetNucName; // add [ for identification as ion /* if( aNucleus ) { G4cout << aNucleus << " NUCLEUS Z " << aNucleus->GetZ_asInt() << " A Z " << aNucleus->GetA_asInt() << G4endl; }*/ secoParticles[targetNucName] = 1; for( ite = secoList->begin(); ite != secoList->end(); ite++ ){ // G4cout << " SECOLIST " << (*ite)->GetDefinition()->GetParticleName() << G4endl; //--- build map of particles G4ParticleDefinition* particle = (*ite)->GetDefinition(); if( bOnlyIons && particle != G4Gamma::Gamma() && particle->GetAtomicNumber() <= 4 ) continue; if( bNotGammas && particle == G4Gamma::Gamma() ) continue; secoParticles[particle->GetParticleName()] ++; } mhmmsi::iterator itemh = theProcessSecoParticles.find( hadProc ); if( itemh == theProcessSecoParticles.end() ) { G4cerr << " PROCESS " << hadProc->GetProcessName(); G4Exception("SHGetHadronicChannelXSUA::UserSteppingAction", "Error", FatalException, "BUG: process not defined in theProcessSecoParticles"); } /* std::map::const_iterator itesp; G4cout << "SECO"; for( itesp = secoParticles.begin(); itesp != secoParticles.end(); itesp++ ) { G4cout << " " << (*itesp).first << " " << (*itesp).second; } G4cout << G4endl; */ // G4cout << " add interaction " << ((*itemh).first)->GetProcessName() << " " << ((*itemh).second).size() << G4endl; ((*itemh).second)[ secoParticles ] ++; } //------------------------------------------------------------------ void SHGetHadronicChannelXSUA::EndOfRunAction(const G4Run*) { mhmmsi::iterator itemh; mmsi::iterator itemp; std::map::reverse_iterator itepi; std::map countPerProcess; std::map< G4String, std::map* > particleYields; // nucleus - (particle - N) //----- Loop to processes for( itemh = theProcessSecoParticles.begin(); itemh != theProcessSecoParticles.end(); itemh++ ) { G4int nEventsProc = 0; std::map::const_iterator iteXS = theCrossSections.find( (*itemh).first ); std::map* >::const_iterator itepys; for( itepys = particleYields.begin(); itepys != particleYields.end(); itepys++ ){ delete (*itepys).second; } particleYields.clear(); //--- Get channel map mmsi secoChannelCount = (*itemh).second; // map of channel - N, where channel is list of particle-N // add channel by ion type std::map countPerChannelAll; // add all channels that are equal without taking into the gammas for( itemp = secoChannelCount.begin(); itemp != secoChannelCount.end(); itemp++ ){ std::map channelParticleList = (*itemp).first; G4cout << "@@SHGetHadronicChannelXSUA " << theCurrentMaterial->GetName() << " " << thePrimaryKineticEnergy << " MeV " << ((*itemh).first)->GetProcessName() << " CHANNEL "; //build channel name G4String channelName = ""; std::vector secoParticles; G4String targetNucName = ""; for(itepi = channelParticleList.rbegin(); itepi != channelParticleList.rend(); itepi++) { G4String partName = (*itepi).first; G4cout << partName << " " << (*itepi).second << " "; if( partName.find("[") == 0 ) targetNucName = partName; if( (*itepi).first != "gamma" ) { for( G4int ip = 0; ip < (*itepi).second; ip++ ) { // if( channelName != "" ) channelName += "_"; // channelName += ((*itepi).first)->GetParticleName(); secoParticles.push_back( partName ); // G4cout << " secoParticles.insert " << (*itepi).first->GetParticleName() << " " << secoParticles.size() << G4endl; } } } if( !thePrimaryParticle ) return; // this run has 0 events (run/beamOn = N* nEventsPerRun) channelName = GmG4Utils::GetInelasticName( secoParticles, thePrimaryParticle ); countPerChannelAll[channelName] += (*itemp).second; countPerProcess[((*itemh).first)] += (*itemp).second; G4cout << " = " << (*itemp).second << G4endl; nEventsProc += (*itemp).second; std::map* >::iterator itepy = particleYields.find(targetNucName); std::map* partYield; if( itepy == particleYields.end() ) { partYield = new std::map(); // G4cout << partYield << " new " << targetNucName << G4endl; particleYields[targetNucName] = partYield; } else { partYield = (*itepy).second; } for(itepi = channelParticleList.rbegin(); itepi != channelParticleList.rend(); itepi++) { G4String partName = (*itepi).first; if( partName != targetNucName ) { (*partYield)[partName] += (*itepi).second * (*itemp).second; // G4cout << targetNucName << " partYield " << partName << " " << (*itepi).second << " " << (*partYield)[partName] << G4endl; } } } std::map::const_iterator iteci; for( iteci = countPerChannelAll.begin(); iteci != countPerChannelAll.end(); iteci++ ){ std::string channelName = (*iteci).first; size_t it = channelName.find("("); std::string targetNucName = channelName.substr(0,it); G4cout << "@@SHGetHadronicChannelXSUA " << theCurrentMaterial->GetName() << " " << thePrimaryKineticEnergy << " MeV" << " " << ((*itemh).first)->GetProcessName() << " CHANNEL_ALL= " << channelName // << " = " << (*iteci).second << "/" << nEventsProc // << "*XS(micro)= " << (*iteci).second*(*iteXS).second/nEventsProc << G4endl; // << " = " << (*iteci).second << "/" << theNEvents[targetNucName] // << "*XS(micro)= " << (*iteci).second*(*iteXS).second/ theNEvents[targetNucName] << G4endl; << " = " << (*iteci).second << "/" << theNEventsInRun << " = " << G4double((*iteci).second)/theNEventsInRun << "*XS(micro)= " << (*iteci).second*(*iteXS).second/ theNEventsInRun << " barn " << G4endl; if( bHistos ) { G4int ichan; std::map::const_iterator itec = theChannels.find(channelName); if( itec == theChannels.end() ) { ichan = theChannels.size(); theChannels[channelName] = ichan; } else { ichan = (*itec).second; } G4int iproc = SHUtils::GetProcessID( (*itemh).first ); std::map::const_iterator item = theMaterials.find(theCurrentMaterial); G4int imat = (*item).second; G4int histoNChann = theHistoNumber + imat*theMateIdxMax*theChanIdxMax*10 + ichan*theChanIdxMax*10 + iproc*10; TH1F* hisNP = theAnaMgr->GetHisto1(histoNChann+3,0); TH1F* hisXSmic = 0; if( !hisNP ) { G4String hname = ((*itemh).first)->GetProcessName() + "_" + theCurrentMaterial->GetName() + "_" + channelName; theAnaMgr->CreateHisto1D((hname+"_NPART").c_str(), theNStepsE, theMinE, theMaxE, histoNChann+3 ); theAnaMgr->CreateHisto1D((hname+"_XSmicro").c_str(), theNStepsE, theMinE, theMaxE, histoNChann+4 ); hisNP = theAnaMgr->GetHisto1(histoNChann+3,0); } hisXSmic = theAnaMgr->GetHisto1(histoNChann+4); hisNP->Fill( thePrimaryKineticEnergy, G4double((*iteci).second)/theNEventsInRun); hisXSmic->Fill( thePrimaryKineticEnergy, (*iteci).second*(*iteXS).second/ theNEventsInRun ); } } if( bParticleYields ) { std::map totalProductYields; // sum for all target nuclei std::map* >::const_iterator itepys; for( itepys = particleYields.begin(); itepys != particleYields.end(); itepys++ ){ G4String targetNucName = (*itepys).first; targetNucName = targetNucName.substr(1,targetNucName.length()); std::map* partYield = (*itepys).second; std::map::const_iterator itepy1; for( itepy1 = partYield->begin(); itepy1 != partYield->end(); itepy1++ ) { G4String productName = (*itepy1).first; G4cout << "@@SHGetHadronicChannelXSUA " << theCurrentMaterial->GetName() << " " << thePrimaryKineticEnergy << " MeV" << " " << ((*itemh).first)->GetProcessName() << " " << targetNucName << " PARTICLE_YIELD " << productName << " = " << (*itepy1).second << "/" << theNEventsInRun << " = " << G4double((*itepy1).second)/theNEventsInRun << " *XS(micro)= " << (*itepy1).second*(*iteXS).second/ theNEventsInRun << G4endl; totalProductYields[productName] += (*itepy1).second; if( bHistos ) { G4int itarg; std::map::const_iterator itec = theTargetNuclei.find(targetNucName); if( itec == theTargetNuclei.end() ) { itarg = theTargetNuclei.size(); theTargetNuclei[targetNucName] = itarg; } else { itarg = (*itec).second; } G4int iproc = SHUtils::GetProcessID( (*itemh).first ); std::map::const_iterator item = theMaterials.find(theCurrentMaterial); G4int imat = (*item).second; G4int iprod; itec = theProducts.find(productName); if( itec == theProducts.end() ) { iprod = theProducts.size(); theProducts[productName] = iprod; } else { iprod = (*itec).second; } G4int histoNYield = theHistoNumber + imat*theMateIdxMax*theChanIdxMax*theProdIdxMax*10 + (itarg+1)*theTargIdxMax*theProdIdxMax*10 + iprod*theProdIdxMax*10 + iproc*10; TH1F* hisY = theAnaMgr->GetHisto1(histoNYield+5,0); TH1F* hisXSmic = 0; if( !hisY ) { G4String hname = ((*itemh).first)->GetProcessName() + "_" + theCurrentMaterial->GetName() + "_" + targetNucName + "_" + productName; theAnaMgr->CreateHisto1D((hname+"_YIELD").c_str(), theNStepsE, theMinE, theMaxE, histoNYield+5 ); theAnaMgr->CreateHisto1D((hname+"_YIELD_XSmicro").c_str(), theNStepsE, theMinE, theMaxE, histoNYield+6 ); hisY = theAnaMgr->GetHisto1(histoNYield+5,0); } hisXSmic = theAnaMgr->GetHisto1(histoNYield+6); hisY->Fill( thePrimaryKineticEnergy, G4double((*itepy1).second)/theNEventsInRun); hisXSmic->Fill( thePrimaryKineticEnergy, (*itepy1).second*(*iteXS).second/ theNEventsInRun ); } } } std::map::const_iterator itepy1; for( itepy1 = totalProductYields.begin(); itepy1 != totalProductYields.end(); itepy1++ ) { G4String productName = (*itepy1).first; G4cout << "@@SHGetHadronicChannelXSUA " << theCurrentMaterial->GetName() << " " << thePrimaryKineticEnergy << " MeV" << " " << ((*itemh).first)->GetProcessName() << " TOTAL_PARTICLE_YIELD " << (*itepy1).first << " = " << (*itepy1).second << "/" << theNEventsInRun << " = " << G4double((*itepy1).second)/theNEventsInRun << " *XS(micro)= " << (*itepy1).second*(*iteXS).second/ theNEventsInRun << G4endl; if( bHistos ) { G4int iproc = SHUtils::GetProcessID( (*itemh).first ); std::map::const_iterator item = theMaterials.find(theCurrentMaterial); G4int imat = (*item).second; std::map::const_iterator itec = theProducts.find(productName); G4int iprod = (*itec).second; G4int histoNYield = theHistoNumber + imat*theMateIdxMax*theChanIdxMax*theProdIdxMax*10 + iprod*theProdIdxMax*10 + iproc*10; TH1F* hisY = theAnaMgr->GetHisto1(histoNYield+5,0); TH1F* hisXSmic = 0; if( !hisY ) { G4String hname = ((*itemh).first)->GetProcessName() + "_" + theCurrentMaterial->GetName() + "_" + productName; theAnaMgr->CreateHisto1D((hname+"_TOTAL_YIELD").c_str(), theNStepsE, theMinE, theMaxE, histoNYield+5 ); theAnaMgr->CreateHisto1D((hname+"_TOTAL_YIELD_XSmicro").c_str(), theNStepsE, theMinE, theMaxE, histoNYield+6 ); hisY = theAnaMgr->GetHisto1(histoNYield+5,0); } hisXSmic = theAnaMgr->GetHisto1(histoNYield+6); hisY->Fill( thePrimaryKineticEnergy, G4double((*itepy1).second)/theNEventsInRun); hisXSmic->Fill( thePrimaryKineticEnergy, (*itepy1).second*(*iteXS).second/ theNEventsInRun ); } } } } /* std::map::const_iterator itecp; for( itecp = countPerProcess.begin(); itecp != countPerProcess.end(); itecp++ ){ std::map::const_iterator iteXS = theCrossSections.find( (*itecp).first ); G4cout << "@@SHGetHadronicChannelXSUA " << theCurrentMaterial->GetName() << " " << thePrimaryKineticEnergy << " PROCESS ALL " << ((*itecp).first)->GetProcessName() << " = " << (*itecp).second << " XS= " << (*iteXS).second << G4endl; } */ theNEvents.clear(); }