Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
CexmcRunManager.cc
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 /*
27  * ============================================================================
28  *
29  * Filename: CexmcRunManager.cc
30  *
31  * Description: run manager
32  *
33  * Version: 1.0
34  * Created: 03.11.2009 20:27:46
35  * Revision: none
36  * Compiler: gcc
37  *
38  * Author: Alexey Radkov (),
39  * Company: PNPI
40  *
41  * ============================================================================
42  */
43 
44 #include <stdlib.h>
45 #include <sys/stat.h>
46 #include <vector>
47 #include <fstream>
48 #ifdef CEXMC_USE_PERSISTENCY
49 #include <boost/archive/binary_iarchive.hpp>
50 #include <boost/archive/binary_oarchive.hpp>
51 #endif
52 #include <G4Eta.hh>
53 #include <G4DigiManager.hh>
54 #include <G4SDManager.hh>
55 #include <G4DecayTable.hh>
56 #include <G4ParticleDefinition.hh>
57 #include <G4ParticleTable.hh>
58 #include <G4UImanager.hh>
59 #include <G4Timer.hh>
60 #include <G4Region.hh>
61 #include <G4RegionStore.hh>
62 #include <G4ProductionCuts.hh>
63 #include <G4VisManager.hh>
64 #include <G4Scene.hh>
65 #include <G4VModel.hh>
66 #include <G4Version.hh>
67 #include "CexmcRunManager.hh"
69 #include "CexmcRunAction.hh"
70 #include "CexmcRun.hh"
71 #include "CexmcPhysicsManager.hh"
72 #include "CexmcProductionModel.hh"
79 #include "CexmcTrackPoints.hh"
83 #include "CexmcEventAction.hh"
84 #include "CexmcParticleGun.hh"
86 #include "CexmcTrackPointsStore.hh"
87 #include "CexmcSetup.hh"
88 #include "CexmcEventSObject.hh"
89 #include "CexmcEventFastSObject.hh"
90 #include "CexmcTrackPointInfo.hh"
91 #include "CexmcEventInfo.hh"
94 #include "CexmcCustomFilterEval.hh"
95 #include "CexmcScenePrimitives.hh"
96 
97 
98 namespace
99 {
100  G4String gdmlFileExtension( ".gdml" );
101  G4String gdmlbz2FileExtension( ".gdml.bz2" );
102 }
103 
104 
106  const G4String & rProject_,
107  G4bool overrideExistingProject ) :
108  basePhysicsUsed( CexmcPMFactoryInstance::GetBasePhysics() ),
109  productionModelType( CexmcUnknownProductionModel ),
110  gdmlFileName( "default.gdml" ), shouldGdmlFileBeValidated( true ),
111  zipGdmlFile( false ), projectsDir( "." ), projectId( projectId_ ),
112  rProject( rProject_ ), guiMacroName( "" ), cfFileName( "" ),
113  eventCountPolicy( CexmcCountAllEvents ),
114  skipInteractionsWithoutEDTonWrite( true ),
115  evDataVerboseLevel( CexmcWriteEventDataOnEveryEDT ),
116  rEvDataVerboseLevel( CexmcWriteNoEventData ), numberOfEventsProcessed( 0 ),
117  numberOfEventsProcessedEffective( 0 ), curEventRead( 0 ),
118 #ifdef CEXMC_USE_PERSISTENCY
119  eventsArchive( NULL ), fastEventsArchive( NULL ),
120 #ifdef CEXMC_USE_CUSTOM_FILTER
121  customFilter( NULL ),
122 #endif
123 #endif
124  physicsManager( NULL ), messenger( NULL )
125 {
126  /* this exception must be caught before creating the object! */
127  if ( rProject != "" && rProject == projectId )
129 
130  const char * projectsDirEnv( getenv( "CEXMC_PROJECTS_DIR" ) );
131 
132  if ( projectsDirEnv )
133  projectsDir = projectsDirEnv;
134 
135  struct stat tmp;
136  if ( ProjectIsSaved() &&
137  stat( ( projectsDir + "/" + projectId + ".rdb" ).c_str(), &tmp ) == 0
138  && ! overrideExistingProject )
140 
141  messenger = new CexmcRunManagerMessenger( this );
142 
143 #ifdef CEXMC_USE_PERSISTENCY
144  if ( ProjectIsRead() )
146 #endif
147 }
148 
149 
151 {
152 #ifdef CEXMC_USE_CUSTOM_FILTER
153  delete customFilter;
154 #endif
155  delete messenger;
156 }
157 
158 
159 #ifdef CEXMC_USE_PERSISTENCY
160 
162 {
163  if ( ! ProjectIsRead() )
164  return;
165 
166  /* read run data */
167  std::ifstream runDataFile( ( projectsDir + "/" + rProject + ".rdb" ).
168  c_str() );
169  if ( ! runDataFile )
171 
172  {
173  boost::archive::binary_iarchive archive( runDataFile );
174  archive >> sObject;
175  }
176 
177  basePhysicsUsed = sObject.basePhysicsUsed;
178 
179  productionModelType = sObject.productionModelType;
180 
181  /* read gdml file */
182  G4String cmd;
183  if ( ProjectIsSaved() )
184  {
185  G4String fileExtension( zipGdmlFile ? gdmlbz2FileExtension :
186  gdmlFileExtension );
187  cmd = G4String( "cp " ) + projectsDir + "/" + rProject + fileExtension +
188  " " + projectsDir + "/" + projectId + fileExtension;
189  if ( system( cmd ) != 0 )
191  }
192 
193  if ( zipGdmlFile )
194  {
195  cmd = G4String( "bunzip2 " ) + projectsDir + "/" + rProject +
196  gdmlbz2FileExtension;
197  if ( system( cmd ) != 0 )
199  }
200 
201  gdmlFileName = projectsDir + "/" + rProject + gdmlFileExtension;
202 }
203 
204 
205 void CexmcRunManager::ReadProject( void )
206 {
207  if ( ! ProjectIsRead() )
208  return;
209 
210  if ( ! physicsManager )
212 
214  sObject.angularRanges );
215  G4DecayTable * etaDecayTable( G4Eta::Definition()->GetDecayTable() );
216  for ( CexmcDecayBranchesStore::const_iterator
217  k( sObject.etaDecayTable.GetDecayBranches().begin() );
218  k != sObject.etaDecayTable.GetDecayBranches().end(); ++k )
219  {
220  etaDecayTable->GetDecayChannel( k->first )->SetBR( k->second );
221  }
222 
224  sObject.fermiMotionIsOn, false );
225  eventCountPolicy = sObject.eventCountPolicy;
226 
227  G4Region * region( G4RegionStore::GetInstance()->GetRegion(
229  if ( ! region )
231 
232  region->GetProductionCuts()->SetProductionCuts(
233  sObject.calorimeterRegCuts );
234 
235  const CexmcPrimaryGeneratorAction * primaryGeneratorAction(
236  static_cast< const CexmcPrimaryGeneratorAction * >(
238  CexmcPrimaryGeneratorAction * thePrimaryGeneratorAction(
239  const_cast< CexmcPrimaryGeneratorAction * >(
240  primaryGeneratorAction ) );
241  CexmcParticleGun * particleGun(
242  thePrimaryGeneratorAction->GetParticleGun() );
243  G4ParticleDefinition * particleDefinition(
244  G4ParticleTable::GetParticleTable()->FindParticle(
245  sObject.beamParticle ) );
246  if ( ! particleDefinition )
248 
249  particleGun->SetParticleDefinition( particleDefinition );
250  particleGun->SetOrigPosition( sObject.beamPos, false );
251  particleGun->SetOrigDirection( sObject.beamDir, false );
252  particleGun->SetOrigMomentumAmp( sObject.beamMomentumAmp, false );
253 
254  thePrimaryGeneratorAction->SetFwhmPosX( sObject.beamFwhmPosX, false );
255  thePrimaryGeneratorAction->SetFwhmPosY( sObject.beamFwhmPosY, false );
256  thePrimaryGeneratorAction->SetFwhmDirX( sObject.beamFwhmDirX, false );
257  thePrimaryGeneratorAction->SetFwhmDirY( sObject.beamFwhmDirY, false );
258  thePrimaryGeneratorAction->SetFwhmMomentumAmp(
259  sObject.beamFwhmMomentumAmp, false );
260 
261  G4DigiManager * digiManager( G4DigiManager::GetDMpointer() );
262  CexmcEnergyDepositDigitizer * edDigitizer(
263  static_cast< CexmcEnergyDepositDigitizer * >(
264  digiManager->FindDigitizerModule( CexmcEDDigitizerName ) ) );
265  if ( ! edDigitizer )
267 
268  edDigitizer->SetMonitorThreshold( sObject.monitorEDThreshold, false );
269  edDigitizer->SetVetoCounterLeftThreshold(
270  sObject.vetoCounterEDLeftThreshold, false );
271  edDigitizer->SetVetoCounterRightThreshold(
272  sObject.vetoCounterEDRightThreshold, false );
273  edDigitizer->SetCalorimeterLeftThreshold(
274  sObject.calorimeterEDLeftThreshold, false );
275  edDigitizer->SetCalorimeterRightThreshold(
276  sObject.calorimeterEDRightThreshold, false );
277  edDigitizer->SetCalorimeterTriggerAlgorithm(
278  sObject.calorimeterTriggerAlgorithm, false );
279  edDigitizer->SetOuterCrystalsVetoAlgorithm(
280  sObject.outerCrystalsVetoAlgorithm, false );
281  edDigitizer->SetOuterCrystalsVetoFraction(
282  sObject.outerCrystalsVetoFraction, false );
283  edDigitizer->ApplyFiniteCrystalResolution(
284  sObject.applyFiniteCrystalResolution, false );
285  edDigitizer->SetCrystalResolutionData( sObject.crystalResolutionData );
286 
287  const CexmcEventAction * eventAction(
288  static_cast< const CexmcEventAction * >( userEventAction ) );
289  if ( ! eventAction )
291 
292  CexmcEventAction * theEventAction( const_cast< CexmcEventAction * >(
293  eventAction ) );
294  CexmcChargeExchangeReconstructor * reconstructor(
295  theEventAction->GetReconstructor() );
296  if ( ! reconstructor )
298 
299  reconstructor->SetCalorimeterEntryPointDefinitionAlgorithm(
300  sObject.epDefinitionAlgorithm );
301  reconstructor->SetCalorimeterEntryPointDepthDefinitionAlgorithm(
302  sObject.epDepthDefinitionAlgorithm );
303  reconstructor->SetCrystalSelectionAlgorithm( sObject.csAlgorithm );
304  reconstructor->UseInnerRefCrystal( sObject.useInnerRefCrystal );
305  reconstructor->SetCalorimeterEntryPointDepth( sObject.epDepth );
306  reconstructor->UseTableMass( sObject.useTableMass );
307  reconstructor->UseMassCut( sObject.useMassCut );
308  reconstructor->SetMassCutOPCenter( sObject.mCutOPCenter );
309  reconstructor->SetMassCutNOPCenter( sObject.mCutNOPCenter );
310  reconstructor->SetMassCutOPWidth( sObject.mCutOPWidth );
311  reconstructor->SetMassCutNOPWidth( sObject.mCutNOPWidth );
312  reconstructor->SetMassCutEllipseAngle( sObject.mCutAngle );
313  reconstructor->UseAbsorbedEnergyCut( sObject.useAbsorbedEnergyCut );
314  reconstructor->SetAbsorbedEnergyCutCLCenter( sObject.aeCutCLCenter );
315  reconstructor->SetAbsorbedEnergyCutCRCenter( sObject.aeCutCRCenter );
316  reconstructor->SetAbsorbedEnergyCutCLWidth( sObject.aeCutCLWidth );
317  reconstructor->SetAbsorbedEnergyCutCRWidth( sObject.aeCutCRWidth );
318  reconstructor->SetAbsorbedEnergyCutEllipseAngle( sObject.aeCutAngle );
319  reconstructor->SetExpectedMomentumAmp( sObject.expectedMomentumAmp );
320  reconstructor->SetEDCollectionAlgorithm( sObject.edCollectionAlgorithm );
321 
322  physicsManager->SetProposedMaxIL( sObject.proposedMaxIL );
323 
324  rEvDataVerboseLevel = sObject.evDataVerboseLevel;
326 }
327 
328 
329 void CexmcRunManager::SaveProject( void )
330 {
331  if ( ! ProjectIsSaved() )
332  return;
333 
334  /* save run data */
335  if ( ! physicsManager )
337 
338  CexmcSimpleDecayTableStore etaDecayTable(
339  G4Eta::Definition()->GetDecayTable() );
340  const CexmcPrimaryGeneratorAction * primaryGeneratorAction(
341  static_cast< const CexmcPrimaryGeneratorAction * >(
343  CexmcPrimaryGeneratorAction * thePrimaryGeneratorAction(
344  const_cast< CexmcPrimaryGeneratorAction * >(
345  primaryGeneratorAction ) );
346  CexmcParticleGun * particleGun(
347  thePrimaryGeneratorAction->GetParticleGun() );
348 
349  G4DigiManager * digiManager( G4DigiManager::GetDMpointer() );
350  CexmcEnergyDepositDigitizer * edDigitizer(
351  static_cast< CexmcEnergyDepositDigitizer * >(
352  digiManager->FindDigitizerModule( CexmcEDDigitizerName ) ) );
353  if ( ! edDigitizer )
355 
356  const CexmcEventAction * eventAction(
357  static_cast< const CexmcEventAction * >( userEventAction ) );
358  CexmcEventAction * theEventAction( const_cast< CexmcEventAction * >(
359  eventAction ) );
360  CexmcChargeExchangeReconstructor * reconstructor(
361  theEventAction->GetReconstructor() );
362 
363  std::vector< G4double > calorimeterRegCuts( 4 );
364  if ( ProjectIsRead() )
365  {
366  calorimeterRegCuts = sObject.calorimeterRegCuts;
367  }
368  else
369  {
370  G4Region * region( G4RegionStore::GetInstance()->GetRegion(
372  if ( ! region )
374 
375  calorimeterRegCuts = region->GetProductionCuts()->GetProductionCuts();
376  }
377 
378  CexmcNmbOfHitsInRanges nmbOfHitsSampled;
379  CexmcNmbOfHitsInRanges nmbOfHitsSampledFull;
380  CexmcNmbOfHitsInRanges nmbOfHitsTriggeredRealRange;
381  CexmcNmbOfHitsInRanges nmbOfHitsTriggeredRecRange;
382  CexmcNmbOfHitsInRanges nmbOfOrphanHits;
383  G4int nmbOfFalseHitsTriggeredEDT( 0 );
384  G4int nmbOfFalseHitsTriggeredRec( 0 );
385  G4int nmbOfSavedEvents( 0 );
386  G4int nmbOfSavedFastEvents( 0 );
387  CexmcRun * run( static_cast< CexmcRun * >( currentRun ) );
388 
389  if ( run )
390  {
391  nmbOfHitsSampled = run->GetNmbOfHitsSampled();
392  nmbOfHitsSampledFull = run->GetNmbOfHitsSampledFull();
393  nmbOfHitsTriggeredRealRange = run->GetNmbOfHitsTriggeredRealRange();
394  nmbOfHitsTriggeredRecRange = run->GetNmbOfHitsTriggeredRecRange();
395  nmbOfOrphanHits = run->GetNmbOfOrphanHits();
396  nmbOfFalseHitsTriggeredEDT = run->GetNmbOfFalseHitsTriggeredEDT();
397  nmbOfFalseHitsTriggeredRec = run->GetNmbOfFalseHitsTriggeredRec();
398  nmbOfSavedEvents = run->GetNmbOfSavedEvents();
399  nmbOfSavedFastEvents = run->GetNmbOfSavedFastEvents();
400  }
401 
402  CexmcRunSObject sObjectToWrite = {
406  calorimeterRegCuts, eventCountPolicy,
407  particleGun->GetParticleDefinition()->GetParticleName(),
408  particleGun->GetOrigPosition(), particleGun->GetOrigDirection(),
409  particleGun->GetOrigMomentumAmp(),
410  primaryGeneratorAction->GetFwhmPosX(),
411  primaryGeneratorAction->GetFwhmPosY(),
412  primaryGeneratorAction->GetFwhmDirX(),
413  primaryGeneratorAction->GetFwhmDirY(),
414  primaryGeneratorAction->GetFwhmMomentumAmp(),
415  edDigitizer->GetMonitorThreshold(),
416  edDigitizer->GetVetoCounterLeftThreshold(),
417  edDigitizer->GetVetoCounterRightThreshold(),
418  edDigitizer->GetCalorimeterLeftThreshold(),
419  edDigitizer->GetCalorimeterRightThreshold(),
420  edDigitizer->GetCalorimeterTriggerAlgorithm(),
421  edDigitizer->GetOuterCrystalsVetoAlgorithm(),
422  edDigitizer->GetOuterCrystalsVetoFraction(),
423  edDigitizer->IsFiniteCrystalResolutionApplied(),
424  edDigitizer->GetCrystalResolutionData(),
425  reconstructor->GetCalorimeterEntryPointDefinitionAlgorithm(),
426  reconstructor->GetCalorimeterEntryPointDepthDefinitionAlgorithm(),
427  reconstructor->GetCrystalSelectionAlgorithm(),
428  reconstructor->IsInnerRefCrystalUsed(),
429  reconstructor->GetCalorimeterEntryPointDepth(),
430  reconstructor->IsTableMassUsed(), reconstructor->IsMassCutUsed(),
431  reconstructor->GetMassCutOPCenter(),
432  reconstructor->GetMassCutNOPCenter(),
433  reconstructor->GetMassCutOPWidth(), reconstructor->GetMassCutNOPWidth(),
434  reconstructor->GetMassCutEllipseAngle(),
435  reconstructor->IsAbsorbedEnergyCutUsed(),
436  reconstructor->GetAbsorbedEnergyCutCLCenter(),
437  reconstructor->GetAbsorbedEnergyCutCRCenter(),
438  reconstructor->GetAbsorbedEnergyCutCLWidth(),
439  reconstructor->GetAbsorbedEnergyCutCRWidth(),
440  reconstructor->GetAbsorbedEnergyCutEllipseAngle(),
441  nmbOfHitsSampled, nmbOfHitsSampledFull, nmbOfHitsTriggeredRealRange,
442  nmbOfHitsTriggeredRecRange, nmbOfOrphanHits, nmbOfFalseHitsTriggeredEDT,
443  nmbOfFalseHitsTriggeredRec, nmbOfSavedEvents, nmbOfSavedFastEvents,
447  reconstructor->GetExpectedMomentumAmp(),
448  reconstructor->GetEDCollectionAlgorithm(), 0 };
449 
450  std::ofstream runDataFile( ( projectsDir + "/" + projectId + ".rdb" ).
451  c_str() );
452 
453  {
454  boost::archive::binary_oarchive archive( runDataFile );
455  archive << sObjectToWrite;
456  }
457 }
458 
459 #endif
460 
461 
463  G4int nSelect )
464 {
465  G4int iEvent( 0 );
466  G4int iEventEffective( 0 );
467 
468  for ( iEvent = 0; iEventEffective < nEvent; ++iEvent )
469  {
470  currentEvent = GenerateEvent( iEvent );
472  CexmcEventInfo * eventInfo( static_cast< CexmcEventInfo * >(
474  switch ( eventCountPolicy )
475  {
476  case CexmcCountAllEvents :
477  ++iEventEffective;
478  break;
480  if ( eventInfo->TpTriggerIsOk() )
481  ++iEventEffective;
482  break;
484  if ( eventInfo->EdTriggerIsOk() )
485  ++iEventEffective;
486  break;
487  default :
488  ++iEventEffective;
489  break;
490  }
492  UpdateScoring();
493  if ( iEvent < nSelect )
496  currentEvent = 0;
497  if ( runAborted )
498  break;
499  }
500 
501  numberOfEventsProcessed = iEvent;
502  numberOfEventsProcessedEffective = iEventEffective;
503 }
504 
505 
506 #ifdef CEXMC_USE_PERSISTENCY
507 
508 void CexmcRunManager::DoReadEventLoop( G4int nEvent )
509 {
510  G4int iEvent( 0 );
511  G4int iEventEffective( 0 );
512  G4int nEventCount( 0 );
513 
514  if ( ! ProjectIsRead() )
515  return;
516 
517  if ( ! physicsManager )
519 
520  CexmcProductionModel * productionModel(
522  if ( ! productionModel )
524 
525  CexmcSetup * setup( static_cast< CexmcSetup * >( userDetector ) );
526  if ( ! setup )
528 
529  CexmcEventSObject evSObject;
530  CexmcEventFastSObject evFastSObject;
531 
532  /* read events data */
533  std::ifstream eventsDataFile(
534  ( projectsDir + "/" + rProject + ".edb" ).c_str() );
535  if ( ! eventsDataFile )
537 
538  boost::archive::binary_iarchive evArchive( eventsDataFile );
539 
540  std::ifstream eventsFastDataFile(
541  ( projectsDir + "/" + rProject + ".fdb" ).c_str() );
542  if ( ! eventsFastDataFile )
544 
545  boost::archive::binary_iarchive evFastArchive( eventsFastDataFile );
546 
547  G4Event event;
548  currentEvent = &event;
549  G4SDManager * sdManager( G4SDManager::GetSDMpointer() );
550  event.SetHCofThisEvent( sdManager->PrepareNewEvent() );
551  G4HCofThisEvent * hcOfThisEvent( event.GetHCofThisEvent() );
552 
553  G4DigiManager * digiManager( G4DigiManager::GetDMpointer() );
554 
555  G4int hcId( digiManager->GetHitsCollectionID(
558  CexmcEnergyDepositCollection * monitorED(
560  hcOfThisEvent->AddHitsCollection( hcId, monitorED );
561  hcId = digiManager->GetHitsCollectionID(
563  "/" + CexmcDetectorTypeName[ CexmcEDDetector ] );
564  CexmcEnergyDepositCollection * vetoCounterED(
566  hcOfThisEvent->AddHitsCollection( hcId, vetoCounterED );
567  hcId = digiManager->GetHitsCollectionID(
569  "/" + CexmcDetectorTypeName[ CexmcEDDetector ] );
570  CexmcEnergyDepositCollection * calorimeterED(
572  hcOfThisEvent->AddHitsCollection( hcId, calorimeterED );
573  hcId = digiManager->GetHitsCollectionID(
577  hcOfThisEvent->AddHitsCollection( hcId, monitorTP );
578  hcId = digiManager->GetHitsCollectionID(
580  "/" + CexmcDetectorTypeName[ CexmcTPDetector ] );
581  CexmcTrackPointsCollection * vetoCounterTP(
583  hcOfThisEvent->AddHitsCollection( hcId, vetoCounterTP );
584  hcId = digiManager->GetHitsCollectionID(
586  "/" + CexmcDetectorTypeName[ CexmcTPDetector ] );
587  CexmcTrackPointsCollection * calorimeterTP(
589  hcOfThisEvent->AddHitsCollection( hcId, calorimeterTP );
590  hcId = digiManager->GetHitsCollectionID(
592  "/" + CexmcDetectorTypeName[ CexmcTPDetector ] );
594  hcOfThisEvent->AddHitsCollection( hcId, targetTP );
595 
596 #ifdef CEXMC_USE_CUSTOM_FILTER
597  if ( customFilter )
598  customFilter->SetAddressedData( &evFastSObject, &evSObject );
599 #endif
600 
601  G4int nmbOfSavedEvents( rEvDataVerboseLevel == CexmcWriteNoEventData ? 0 :
602  sObject.nmbOfSavedFastEvents );
603  G4bool eventDataWrittenOnEveryTPT( rEvDataVerboseLevel ==
605 
606  for ( int i( 0 ); i < nmbOfSavedEvents; ++i )
607  {
608  evFastArchive >> evFastSObject;
609 
610  if ( nEventCount < curEventRead )
611  {
612  if ( eventDataWrittenOnEveryTPT ||
613  evFastSObject.edDigitizerHasTriggered )
614  {
615  evArchive >> evSObject;
616  if ( evFastSObject.edDigitizerHasTriggered )
617  ++nEventCount;
618  }
619  continue;
620  }
621 
622  ++iEvent;
623 
624  productionModel->SetTriggeredAngularRanges(
625  evFastSObject.opCosThetaSCM );
626  const CexmcAngularRangeList & triggeredAngularRanges(
627  productionModel->GetTriggeredAngularRanges() );
628 
629  if ( ! eventDataWrittenOnEveryTPT &&
630  ! evFastSObject.edDigitizerHasTriggered )
631  {
632 #ifdef CEXMC_USE_CUSTOM_FILTER
633  /* user must be aware that using tpt commands in custom filter
634  * scripts for poor event data sets can easily lead to logical
635  * errors! This is because most of tpt data is only available for
636  * events with EDT trigger. There is no such problem if event data
637  * was written on every TPT event. */
638  if ( customFilter && ! customFilter->EvalTPT() )
639  continue;
640 #endif
641  SaveCurrentTPTEvent( evFastSObject, triggeredAngularRanges,
642  ProjectIsSaved() && ! skipInteractionsWithoutEDTonWrite );
643  continue;
644  }
645 
646  evArchive >> evSObject;
647 
648  G4bool skipEDTOnThisEvent( false );
649 
650 #ifdef CEXMC_USE_CUSTOM_FILTER
651  if ( customFilter && ! customFilter->EvalTPT() )
652  continue;
653  if ( customFilter && ! customFilter->EvalEDT() )
654  {
655  if ( ! eventDataWrittenOnEveryTPT )
656  {
657  SaveCurrentTPTEvent( evFastSObject, triggeredAngularRanges,
658  ProjectIsSaved() );
659  continue;
660  }
661  skipEDTOnThisEvent = true;
662  }
663 #endif
664 
665  event.SetEventID( evSObject.eventId );
666 
667  /* AAA: beware! If something throws an exception between AAA and CCC
668  * then there would be segmentation violation in ~THitsMap() as it
669  * would try to delete local variable evSObject's fields.
670  * See also BBB */
671 
672  monitorED->GetMap()->operator[]( 0 ) = &evSObject.monitorED;
673  vetoCounterED->GetMap()->operator[]( 0 ) = &evSObject.vetoCounterEDLeft;
674  vetoCounterED->GetMap()->operator[]( 1 <<
676  &evSObject.vetoCounterEDRight;
677  G4int row( 0 );
678  G4int column( 0 );
679  for ( CexmcEnergyDepositCalorimeterCollection::iterator
680  k( evSObject.calorimeterEDLeftCollection.begin() );
681  k != evSObject.calorimeterEDLeftCollection.end(); ++k )
682  {
683  G4int index( row <<
685  column = 0;
686  for ( CexmcEnergyDepositCrystalRowCollection::iterator
687  l( k->begin() ); l != k->end(); ++l )
688  {
689  calorimeterED->GetMap()->operator[]( index | column ) = &*l;
690  ++column;
691  }
692  ++row;
693  }
694  row = 0;
695  for ( CexmcEnergyDepositCalorimeterCollection::iterator
696  k( evSObject.calorimeterEDRightCollection.begin() );
697  k != evSObject.calorimeterEDRightCollection.end(); ++k )
698  {
699  G4int index(
701  | row <<
703  column = 0;
704  for ( CexmcEnergyDepositCrystalRowCollection::iterator
705  l( k->begin() ); l != k->end(); ++l )
706  {
707  calorimeterED->GetMap()->operator[]( index | column ) = &*l;
708  ++column;
709  }
710  ++row;
711  }
712 
713  CexmcTrackPointInfo monitorTPInfo( evSObject.monitorTP );
714  CexmcTrackPointInfo targetTPBeamParticleInfo(
715  evSObject.targetTPBeamParticle );
716  CexmcTrackPointInfo targetTPOutputParticleInfo(
717  evSObject.targetTPOutputParticle );
718  CexmcTrackPointInfo targetTPNucleusParticleInfo(
719  evSObject.targetTPNucleusParticle );
720  CexmcTrackPointInfo targetTPOutputParticleDecayProductParticle1Info(
721  evSObject.targetTPOutputParticleDecayProductParticle1 );
722  CexmcTrackPointInfo targetTPOutputParticleDecayProductParticle2Info(
723  evSObject.targetTPOutputParticleDecayProductParticle2 );
724  CexmcTrackPointInfo vetoCounterTPLeftInfo(
725  evSObject.vetoCounterTPLeft );
726  CexmcTrackPointInfo vetoCounterTPRightInfo(
727  evSObject.vetoCounterTPRight );
728  CexmcTrackPointInfo calorimeterTPLeftInfo(
729  evSObject.calorimeterTPLeft );
730  CexmcTrackPointInfo calorimeterTPRightInfo(
731  evSObject.calorimeterTPRight );
732 
733  if ( monitorTPInfo.IsValid() )
734  monitorTP->GetMap()->operator[]( monitorTPInfo.trackId ) =
735  &monitorTPInfo;
736  if ( targetTPBeamParticleInfo.IsValid() )
737  targetTP->GetMap()->operator[](
738  targetTPBeamParticleInfo.trackId ) =
739  &targetTPBeamParticleInfo;
740  if ( targetTPOutputParticleInfo.IsValid() )
741  targetTP->GetMap()->operator[](
742  targetTPOutputParticleInfo.trackId ) =
743  &targetTPOutputParticleInfo;
744  if ( targetTPNucleusParticleInfo.IsValid() )
745  targetTP->GetMap()->operator[](
746  targetTPNucleusParticleInfo.trackId ) =
747  &targetTPNucleusParticleInfo;
748  if ( targetTPOutputParticleDecayProductParticle1Info.IsValid() )
749  targetTP->GetMap()->operator[](
750  targetTPOutputParticleDecayProductParticle1Info.trackId ) =
751  &targetTPOutputParticleDecayProductParticle1Info;
752  if ( targetTPOutputParticleDecayProductParticle2Info.IsValid() )
753  targetTP->GetMap()->operator[](
754  targetTPOutputParticleDecayProductParticle2Info.trackId ) =
755  &targetTPOutputParticleDecayProductParticle2Info;
756  if ( vetoCounterTPLeftInfo.IsValid() )
757  vetoCounterTP->GetMap()->operator[](
758  vetoCounterTPLeftInfo.trackId ) = &vetoCounterTPLeftInfo;
759  if ( vetoCounterTPRightInfo.IsValid() )
760  vetoCounterTP->GetMap()->operator[](
762  vetoCounterTPRightInfo.trackId ) = &vetoCounterTPRightInfo;
763 
765  if ( calorimeterTPLeftInfo.IsValid() )
766  {
767  pos = calorimeterTPLeftInfo.positionLocal;
768  setup->ConvertToCrystalGeometry(
769  calorimeterTPLeftInfo.positionLocal, row, column, pos );
770  calorimeterTPLeftInfo.positionLocal = pos;
771  calorimeterTP->GetMap()->operator[](
776  calorimeterTPLeftInfo.trackId ) = &calorimeterTPLeftInfo;
777  }
778  if ( calorimeterTPRightInfo.IsValid() )
779  {
780  pos = calorimeterTPRightInfo.positionLocal;
781  setup->ConvertToCrystalGeometry(
782  calorimeterTPRightInfo.positionLocal, row, column, pos );
783  calorimeterTPRightInfo.positionLocal = pos;
784  calorimeterTP->GetMap()->operator[](
790  calorimeterTPRightInfo.trackId ) = &calorimeterTPRightInfo;
791  }
792 
793  productionModel->SetProductionModelData(
794  evSObject.productionModelData );
795 
796  const CexmcEventAction * eventAction(
797  static_cast< const CexmcEventAction * >( userEventAction ) );
798  if ( ! eventAction )
799  {
800  /* BBB: all hits collections must be cleared before throwing
801  * anything from here, otherwise ~THitsMap() will try to delete
802  * local variable evSObject's fields like monitorED etc. */
803  monitorED->GetMap()->clear();
804  vetoCounterED->GetMap()->clear();
805  calorimeterED->GetMap()->clear();
806  monitorTP->GetMap()->clear();
807  targetTP->GetMap()->clear();
808  vetoCounterTP->GetMap()->clear();
809  calorimeterTP->GetMap()->clear();
811  }
812 
813  if ( skipEDTOnThisEvent )
814  event.SetUserInformation( new CexmcEventInfo( false, false,
815  false ) );
816 
817  CexmcEventAction * theEventAction( const_cast< CexmcEventAction * >(
818  eventAction ) );
819  theEventAction->EndOfEventAction( &event );
820 
821  CexmcEventInfo * eventInfo( static_cast< CexmcEventInfo * >(
822  event.GetUserInformation() ) );
823 
824  if ( eventInfo->EdTriggerIsOk() )
825  ++iEventEffective;
826 
827  delete eventInfo;
828  event.SetUserInformation( NULL );
829 
830  monitorED->GetMap()->clear();
831  vetoCounterED->GetMap()->clear();
832  calorimeterED->GetMap()->clear();
833  monitorTP->GetMap()->clear();
834  targetTP->GetMap()->clear();
835  vetoCounterTP->GetMap()->clear();
836  calorimeterTP->GetMap()->clear();
837 
838  /* CCC: see AAA */
839 
840  if ( nEvent > 0 && iEventEffective == nEvent )
841  break;
842  }
843 
844  curEventRead = nEventCount + iEventEffective;
845 
846  numberOfEventsProcessed = iEvent;
847  numberOfEventsProcessedEffective = iEventEffective;
848 
849 #ifdef CEXMC_USE_CUSTOM_FILTER
850  if ( customFilter )
851  customFilter->SetAddressedData( NULL, NULL );
852 #endif
853 }
854 
855 
856 void CexmcRunManager::SaveCurrentTPTEvent(
857  const CexmcEventFastSObject & evFastSObject,
858  const CexmcAngularRangeList & angularRanges,
859  G4bool writeToDatabase )
860 {
861  CexmcRun * run( static_cast< CexmcRun * >( currentRun ) );
862 
863  if ( ! run )
864  return;
865 
866  for ( CexmcAngularRangeList::const_iterator k( angularRanges.begin() );
867  k != angularRanges.end(); ++k )
868  {
869  run->IncrementNmbOfHitsSampledFull( k->index );
870  if ( evFastSObject.edDigitizerMonitorHasTriggered )
871  run->IncrementNmbOfHitsSampled( k->index );
872  }
873 
874  if ( writeToDatabase )
875  {
876  fastEventsArchive->operator<<( evFastSObject );
877  run->IncrementNmbOfSavedFastEvents();
878  }
879 }
880 
881 #endif
882 
883 
884 /* mostly adopted from G4RunManager::DoEventLoop() */
885 void CexmcRunManager::DoEventLoop( G4int nEvent, const char * macroFile,
886  G4int nSelect )
887 {
888  if ( verboseLevel > 0 )
889  timer->Start();
890 
891  G4String cmd;
892  if ( macroFile != 0 )
893  {
894  if ( nSelect < 0 )
895  nSelect = nEvent;
896  cmd = "/control/execute ";
897  cmd += macroFile;
898  }
899  else
900  {
901  nSelect = -1;
902  }
903 
904  numberOfEventsProcessed = 0;
905  numberOfEventsProcessedEffective = 0;
906 
907 #ifdef CEXMC_USE_PERSISTENCY
908  eventsArchive = NULL;
909  fastEventsArchive = NULL;
910  if ( ProjectIsRead() )
911  {
912  if ( ProjectIsSaved() )
913  {
914  std::ofstream eventsDataFile(
915  ( projectsDir + "/" + projectId + ".edb" ).c_str() );
916  boost::archive::binary_oarchive eventsArchive_( eventsDataFile );
917  std::ofstream fastEventsDataFile(
918  ( projectsDir + "/" + projectId + ".fdb" ).c_str() );
919  boost::archive::binary_oarchive fastEventsArchive_(
920  fastEventsDataFile );
921  eventsArchive = &eventsArchive_;
922  fastEventsArchive = &fastEventsArchive_;
923  DoReadEventLoop( nEvent );
924  }
925  else
926  {
927  DoReadEventLoop( nEvent );
928  }
929  }
930  else
931  {
932  if ( ProjectIsSaved() )
933  {
934  std::ofstream eventsDataFile(
935  ( projectsDir + "/" + projectId + ".edb" ).c_str() );
936  boost::archive::binary_oarchive eventsArchive_( eventsDataFile );
937  std::ofstream fastEventsDataFile(
938  ( projectsDir + "/" + projectId + ".fdb" ).c_str() );
939  boost::archive::binary_oarchive fastEventsArchive_(
940  fastEventsDataFile );
941  eventsArchive = &eventsArchive_;
942  fastEventsArchive = &fastEventsArchive_;
943  DoCommonEventLoop( nEvent, cmd, nSelect );
944  }
945  else
946  {
947  DoCommonEventLoop( nEvent, cmd, nSelect );
948  }
949  }
950  eventsArchive = NULL;
951  fastEventsArchive = NULL;
952 #else
953  DoCommonEventLoop( nEvent, cmd, nSelect );
954 #endif
955 
956  if ( verboseLevel > 0 )
957  {
958  timer->Stop();
959  G4cout << "Run terminated." << G4endl;
960  G4cout << "Run Summary" << G4endl;
961  if ( runAborted )
962  {
963  G4cout << " Run Aborted after " << numberOfEventsProcessed <<
964  " events processed." << G4endl;
965  }
966  else
967  {
968  G4cout << " Number of events processed : " <<
969  numberOfEventsProcessed << ", effectively: " <<
970  numberOfEventsProcessedEffective << G4endl;
971  }
972  G4cout << " " << *timer << G4endl;
973  }
974 }
975 
976 
977 #ifdef CEXMC_USE_PERSISTENCY
978 
979 void CexmcRunManager::PrintReadRunData( void ) const
980 {
981  if ( ! ProjectIsRead() )
982  return;
983 
984  G4bool refCrystalInfoPrinted( false );
985 
986  G4cout << CEXMC_LINE_START << "Run data read from project '" << rProject <<
987  "'" << G4endl;
988  G4cout << " (archive class version " <<
989  sObject.actualVersion << ")" << G4endl;
990  if ( ! sObject.rProject.empty() )
991  {
992  G4cout << " -- Based on project '" << sObject.rProject << "'" <<
993  G4endl;
994  if ( ! sObject.cfFileName.empty() )
995  G4cout << " -- Custom filter script '" << sObject.cfFileName <<
996  "' was used" << G4endl;
997  }
998  G4cout << " -- Event data verbose level (0 - not saved, 1 - triggers, "
999  "2 - interactions): " << sObject.evDataVerboseLevel << G4endl;
1000  if ( ! sObject.rProject.empty() )
1001  {
1002  if ( sObject.evDataVerboseLevel == CexmcWriteEventDataOnEveryEDT )
1003  {
1004  G4cout << " -- (fdb file contains " <<
1005  ( sObject.interactionsWithoutEDTWereSkipped ?
1006  "only interactions when an event was triggered" :
1007  "all interactions" ) << ")" << std::endl;
1008  }
1009  }
1010  G4cout << " -- Base physics used"
1011  "(1 - QGSP_BERT, 2 - QGSP_BIC_EMY, 3 - FTFP_BERT): " <<
1012  sObject.basePhysicsUsed << G4endl;
1013  G4cout << " -- Production model (1 - pi0, 2 - eta): " <<
1014  sObject.productionModelType << G4endl;
1015  G4cout << " -- Geometry definition file: " << sObject.gdmlFileName <<
1016  G4endl;
1017  G4cout << " -- Angular ranges: " << sObject.angularRanges << G4endl;
1018  G4cout << " -- Eta decay modes: " << G4endl;
1020  G4cout << " -- Fermi motion status (0 - disabled, 1 - enabled): " <<
1021  sObject.fermiMotionIsOn << G4endl;
1022  if ( sObject.calorimeterRegCuts.size() < 4 )
1024  G4cout << " -- Production cuts in calorimeter (gamma, e-, e+, p): " <<
1025  G4BestUnit( sObject.calorimeterRegCuts[ 0 ], "Length" ) << ", " <<
1026  G4BestUnit( sObject.calorimeterRegCuts[ 1 ], "Length" ) << ", " <<
1027  G4BestUnit( sObject.calorimeterRegCuts[ 2 ], "Length" ) << ", " <<
1028  G4BestUnit( sObject.calorimeterRegCuts[ 3 ], "Length" ) << G4endl;
1029  G4cout << " -- Proposed max interaction length in the target: " <<
1030  G4BestUnit( sObject.proposedMaxIL, "Length" ) << G4endl;
1031  G4cout << " -- Event count policy (0 - all, 1 - interaction, 2 - trigger)"
1032  ": " << sObject.eventCountPolicy << G4endl;
1033  G4cout << " -- Number of events (processed / effective / ordered): " <<
1034  sObject.numberOfEventsProcessed << " / " <<
1035  sObject.numberOfEventsProcessedEffective << " / " <<
1036  sObject.numberOfEventsToBeProcessed << G4endl;
1037  G4cout << " -- Incident beam particle: " << sObject.beamParticle << G4endl;
1038  G4cout << " position: " <<
1039  G4BestUnit( sObject.beamPos, "Length" ) << G4endl;
1040  G4cout << " direction: " <<
1041  G4ThreeVector( sObject.beamDir ) << G4endl;
1042  G4cout << " momentum: " <<
1043  G4BestUnit( sObject.beamMomentumAmp, "Energy" ) << G4endl;
1044  G4cout << " momentum fwhm: " << sObject.beamFwhmMomentumAmp <<
1045  G4endl;
1046  G4cout << " pos fwhm (x): " <<
1047  G4BestUnit( sObject.beamFwhmPosX, "Length" ) << G4endl;
1048  G4cout << " pos fwhm (y): " <<
1049  G4BestUnit( sObject.beamFwhmPosY, "Length" ) << G4endl;
1050  G4cout << " dir fwhm (x): " << sObject.beamFwhmDirX / deg <<
1051  " deg" << G4endl;
1052  G4cout << " dir fwhm (y): " << sObject.beamFwhmDirY / deg <<
1053  " deg" << G4endl;
1054  G4cout << " -- Monitor ED threshold: " <<
1055  G4BestUnit( sObject.monitorEDThreshold, "Energy" ) << G4endl;
1056  G4cout << " -- Veto counter (l/r) ED threshold: " <<
1057  G4BestUnit( sObject.vetoCounterEDLeftThreshold, "Energy" ) <<
1058  " / " <<
1059  G4BestUnit( sObject.vetoCounterEDRightThreshold, "Energy" ) <<
1060  G4endl;
1061  G4cout << " -- Calorimeter (l/r) ED threshold: " <<
1062  G4BestUnit( sObject.calorimeterEDLeftThreshold, "Energy" ) <<
1063  " / " <<
1064  G4BestUnit( sObject.calorimeterEDRightThreshold, "Energy" ) <<
1065  G4endl;
1066  G4cout << " -- Calorimeter trigger algorithm (0 - all, 1 - inner): " <<
1067  sObject.calorimeterTriggerAlgorithm << G4endl;
1068  G4cout << " -- Outer crystals veto algorithm "
1069  "(0 - none, 1 - max, 2 - fraction): " <<
1070  sObject.outerCrystalsVetoAlgorithm << G4endl;
1071  if ( sObject.outerCrystalsVetoAlgorithm ==
1073  {
1074  G4cout << " -- Outer crystals veto fraction: " <<
1075  sObject.outerCrystalsVetoFraction << G4endl;
1076  }
1077  G4cout << " -- Finite crystal resolution applied (0 - no, 1 - yes): " <<
1078  sObject.applyFiniteCrystalResolution << G4endl;
1079  if ( sObject.applyFiniteCrystalResolution )
1080  {
1081  G4cout << " -- Crystal resolution data: " <<
1082  sObject.crystalResolutionData;
1083  }
1084  G4cout << " -- Reconstructor settings: " << G4endl;
1085  if ( sObject.expectedMomentumAmp > 0 )
1086  {
1087  G4cout << " -- expected momentum in the target: " <<
1088  G4BestUnit( sObject.expectedMomentumAmp, "Energy" ) << G4endl;
1089  }
1090  G4cout << " -- ed collection algorithm (0 - all, 1 - adjacent): " <<
1091  sObject.edCollectionAlgorithm << G4endl;
1092  if ( sObject.edCollectionAlgorithm == CexmcCollectEDInAdjacentCrystals )
1093  {
1094  G4cout <<
1095  " -- inner crystal used as reference (0 - no, 1 - yes): " <<
1096  sObject.useInnerRefCrystal << G4endl;
1097  refCrystalInfoPrinted = true;
1098  }
1099  G4cout << " -- entry point definition algorithm " << G4endl;
1100  G4cout << " (0 - center of calorimeter, 1 - center of crystal with "
1101  "max ED," << G4endl;
1102  G4cout << " 2 - linear, 3 - square): " <<
1103  sObject.epDefinitionAlgorithm << G4endl;
1104  G4cout << " -- entry point depth definition algorithm "
1105  "(0 - plain, 1 - sphere): " <<
1106  sObject.epDepthDefinitionAlgorithm << G4endl;
1107  G4cout << " -- entry point depth: " <<
1108  G4BestUnit( sObject.epDepth, "Length" ) << G4endl;
1109  if ( sObject.epDefinitionAlgorithm == CexmcEntryPointByLinearEDWeights ||
1110  sObject.epDefinitionAlgorithm == CexmcEntryPointBySqrtEDWeights )
1111  {
1112  G4cout <<
1113  " -- crystal selection algorithm (0 - all, 1 - adjacent): " <<
1114  sObject.csAlgorithm << G4endl;
1115  }
1116  if ( ! refCrystalInfoPrinted &&
1117  ( sObject.epDefinitionAlgorithm ==
1119  ( ( sObject.epDefinitionAlgorithm == CexmcEntryPointBySqrtEDWeights ||
1120  sObject.epDefinitionAlgorithm ==
1122  sObject.csAlgorithm == CexmcSelectAdjacentCrystals ) ) )
1123  {
1124  G4cout <<
1125  " -- inner crystal used as reference (0 - no, 1 - yes): " <<
1126  sObject.useInnerRefCrystal << G4endl;
1127  }
1128  G4cout << " -- table mass of output particle used "
1129  "(0 - no, 1 - yes): " << sObject.useTableMass << G4endl;
1130  G4cout << " -- mass cut is enabled (0 - no, 1 - yes): " <<
1131  sObject.useMassCut << G4endl;
1132  if ( sObject.useMassCut )
1133  {
1134  G4cout << " -- mass cut output particle center: " <<
1135  G4BestUnit( sObject.mCutOPCenter, "Energy" ) << G4endl;
1136  G4cout << " -- mass cut nucleus output particle center: " <<
1137  G4BestUnit( sObject.mCutNOPCenter, "Energy" ) << G4endl;
1138  G4cout << " -- mass cut output particle width of the ellipse: " <<
1139  G4BestUnit( sObject.mCutOPWidth, "Energy" ) << G4endl;
1140  G4cout << " -- mass cut nucleus output particle width of the "
1141  "ellipse: "
1142  << G4BestUnit( sObject.mCutNOPWidth, "Energy" ) << G4endl;
1143  G4cout << " -- mass cut angle of the ellipse: " <<
1144  sObject.mCutAngle / deg << " deg" << G4endl;
1145  }
1146  G4cout << " -- absorbed energy cut is enabled (0 - no, 1 - yes): " <<
1147  sObject.useAbsorbedEnergyCut << G4endl;
1148  if ( sObject.useAbsorbedEnergyCut )
1149  {
1150  G4cout << " -- absorbed energy cut left calorimeter center: " <<
1151  G4BestUnit( sObject.aeCutCLCenter, "Energy" ) << G4endl;
1152  G4cout << " -- absorbed energy cut right calorimeter center: " <<
1153  G4BestUnit( sObject.aeCutCRCenter, "Energy" ) << G4endl;
1154  G4cout << " -- absorbed energy cut left calorimeter width of the "
1155  "ellipse: " <<
1156  G4BestUnit( sObject.aeCutCLWidth, "Energy" ) << G4endl;
1157  G4cout << " -- absorbed energy cut right calorimeter width of the "
1158  "ellipse: "
1159  << G4BestUnit( sObject.aeCutCRWidth, "Energy" ) << G4endl;
1160  G4cout << " -- absorbed energy cut angle of the ellipse: " <<
1161  sObject.aeCutAngle / deg << " deg" << G4endl;
1162  }
1163  G4cout << G4endl;
1164  CexmcRunAction::PrintResults( sObject.nmbOfHitsSampled,
1165  sObject.nmbOfHitsSampledFull,
1166  sObject.nmbOfHitsTriggeredRealRange,
1167  sObject.nmbOfHitsTriggeredRecRange,
1168  sObject.nmbOfOrphanHits,
1169  sObject.angularRanges,
1170  sObject.nmbOfFalseHitsTriggeredEDT,
1171  sObject.nmbOfFalseHitsTriggeredRec );
1172  G4cout << G4endl;
1173 }
1174 
1175 
1176 void CexmcRunManager::ReadAndPrintEventsData( void ) const
1177 {
1178  if ( ! ProjectIsRead() )
1179  return;
1180 
1181  CexmcEventSObject evSObject;
1182 
1183  /* read events data */
1184  std::ifstream eventsDataFile(
1185  ( projectsDir + "/" + rProject + ".edb" ).c_str() );
1186  if ( ! eventsDataFile )
1188 
1189  boost::archive::binary_iarchive evArchive( eventsDataFile );
1190 
1191  for ( int i( 0 ); i < sObject.nmbOfSavedEvents; ++i )
1192  {
1193  evArchive >> evSObject;
1194 
1195  if ( ! evSObject.edDigitizerMonitorHasTriggered )
1196  continue;
1197 
1198  CexmcEnergyDepositStore edStore( evSObject.monitorED,
1199  evSObject.vetoCounterEDLeft, evSObject.vetoCounterEDRight,
1200  evSObject.calorimeterEDLeft, evSObject.calorimeterEDRight,
1201  0, 0, 0, 0, evSObject.calorimeterEDLeftCollection,
1202  evSObject.calorimeterEDRightCollection );
1203 
1204  CexmcTrackPointsStore tpStore( evSObject.monitorTP,
1205  evSObject.targetTPBeamParticle, evSObject.targetTPOutputParticle,
1206  evSObject.targetTPNucleusParticle,
1207  evSObject.targetTPOutputParticleDecayProductParticle1,
1208  evSObject.targetTPOutputParticleDecayProductParticle2,
1209  evSObject.vetoCounterTPLeft, evSObject.vetoCounterTPRight,
1210  evSObject.calorimeterTPLeft, evSObject.calorimeterTPRight );
1211 
1212  const CexmcProductionModelData & pmData(
1213  evSObject.productionModelData );
1214 
1215  G4cout << "Event " << evSObject.eventId << G4endl;
1217  G4cout << " --- Production model data: " << pmData;
1219  }
1220 }
1221 
1222 
1223 void CexmcRunManager::PrintReadData(
1224  const CexmcOutputDataTypeSet & outputData ) const
1225 {
1226  if ( ! ProjectIsRead() )
1227  return;
1228 
1229  G4bool addSpace( false );
1230 
1231  CexmcOutputDataTypeSet::const_iterator found(
1232  outputData.find( CexmcOutputGeometry ) );
1233  if ( found != outputData.end() )
1234  {
1235  G4String cmd( G4String( "cat " ) + projectsDir + "/" + rProject +
1236  gdmlFileExtension );
1237  if ( system( cmd ) != 0 )
1239 
1240  if ( zipGdmlFile )
1241  {
1242  cmd = G4String( "bzip2 " ) + projectsDir + "/" + rProject +
1243  gdmlFileExtension;
1244  if ( system( cmd ) != 0 )
1246  }
1247 
1248  addSpace = true;
1249  }
1250 
1251  found = outputData.find( CexmcOutputEvents );
1252  if ( found != outputData.end() )
1253  {
1254  if ( addSpace )
1255  G4cout << G4endl << G4endl;
1256 
1257  ReadAndPrintEventsData();
1258 
1259  addSpace = true;
1260  }
1261 
1262  found = outputData.find( CexmcOutputRun );
1263  if ( found != outputData.end() )
1264  {
1265  if ( addSpace )
1266  G4cout << G4endl << G4endl;
1267 
1268  G4DecayTable * etaDecayTable( G4Eta::Definition()->GetDecayTable() );
1269  for ( CexmcDecayBranchesStore::const_iterator
1270  k( sObject.etaDecayTable.GetDecayBranches().begin() );
1271  k != sObject.etaDecayTable.GetDecayBranches().end(); ++k )
1272  {
1273  etaDecayTable->GetDecayChannel( k->first )->SetBR( k->second );
1274  }
1275 
1276  PrintReadRunData();
1277  }
1278 }
1279 
1280 
1281 #ifdef CEXMC_USE_CUSTOM_FILTER
1282 
1283 void CexmcRunManager::SetCustomFilter( const G4String & cfFileName_ )
1284 {
1285  if ( customFilter )
1286  {
1287  delete customFilter;
1288  customFilter = NULL;
1289  }
1290 
1291  if ( cfFileName_.empty() )
1292  return;
1293 
1294  /* should not get here */
1295  if ( ! ProjectIsRead() )
1297 
1298  cfFileName = cfFileName_;
1299 
1300  customFilter = new CexmcCustomFilterEval( cfFileName );
1301 }
1302 
1303 #endif
1304 
1305 #endif
1306 
1307 
1309 {
1310  G4VisManager * visManager( static_cast< G4VisManager * >(
1312  if ( ! visManager )
1313  return;
1314 
1315  G4Scene * curScene( visManager->GetCurrentScene() );
1316  if ( ! curScene )
1317  return;
1318 
1319  /* G4Scene declarations lack this kind of typedef */
1320 #if G4VERSION_NUMBER < 960
1321  typedef std::vector< G4VModel * > MList;
1322 #else
1323  typedef std::vector< G4Scene::Model > MList;
1324 #endif
1325  const MList & mList( curScene->GetRunDurationModelList() );
1326 
1327  for ( MList::const_iterator k( mList.begin() ); k != mList.end(); ++k )
1328  {
1329 #if G4VERSION_NUMBER < 960
1330  const G4String & modelDesc( ( *k )->GetGlobalDescription() );
1331 #else
1332  const G4String & modelDesc( k->fpModel->GetGlobalDescription() );
1333 #endif
1334  if ( modelDesc == CexmcScenePrimitivesDescription )
1335  return;
1336  }
1337 
1338  CexmcSetup * setup( static_cast< CexmcSetup * >( userDetector ) );
1339  if ( ! setup )
1341 
1342  /* BEWARE: looks like G4Scene won't delete models from its lists upon
1343  * termination! Hence destructor of the new model won't be called */
1344  curScene->AddRunDurationModel( new CexmcScenePrimitives( setup ) );
1345 }
1346 
1347 
1349 {
1350  const CexmcEventAction * eventAction(
1351  static_cast< const CexmcEventAction * >( userEventAction ) );
1352  if ( ! eventAction )
1354 
1355  CexmcEventAction * theEventAction( const_cast< CexmcEventAction * >(
1356  eventAction ) );
1357  theEventAction->BeamParticleChangeHook();
1358 }
1359 
1360 
1362 {
1363 #ifdef CEXMC_USE_PERSISTENCY
1364  /* save gdml file */
1365  G4String cmd( "" );
1366  CexmcExceptionType exceptionType( CexmcSystemException );
1367 
1368  if ( zipGdmlFile )
1369  {
1370  if ( ProjectIsRead() )
1371  {
1372  cmd = G4String( "bzip2 " ) + projectsDir + "/" + rProject +
1373  gdmlFileExtension;
1374  }
1375  else
1376  {
1377  if ( ProjectIsSaved() )
1378  cmd = G4String( "bzip2 -c " ) + gdmlFileName + " > " +
1379  projectsDir + "/" + projectId + gdmlbz2FileExtension;
1380  }
1381  exceptionType = CexmcFileCompressException;
1382  }
1383  else
1384  {
1385  if ( ! ProjectIsRead() && ProjectIsSaved() )
1386  cmd = G4String( "cp " ) + gdmlFileName + " " + projectsDir + "/" +
1387  projectId + gdmlFileExtension;
1388  /* else already saved in ReadPreinitProjectData() */
1389  }
1390 
1391  if ( ! cmd.empty() && system( cmd ) != 0 )
1392  throw CexmcException( exceptionType );
1393 #endif
1394 }
1395 
static G4DigiManager * GetDMpointer()
void BeamParticleChangeHook(void)
G4VUserDetectorConstruction * userDetector
CexmcExceptionType
void DoEventLoop(G4int nEvent, const char *macroFile, G4int nSelect)
CLHEP::Hep3Vector G4ThreeVector
const G4String CexmcEDDigitizerName("EDDig")
void DoCommonEventLoop(G4int nEvent, const G4String &cmd, G4int nSelect)
G4Event * currentEvent
const G4String CexmcDetectorRoleName[CexmcNumberOfDetectorRoles]
G4DecayTable * GetDecayTable() const
system("rm -rf microbeam.root")
static const G4double pos
static G4ParticleTable * GetParticleTable()
CexmcPhysicsManager * physicsManager
const G4String CexmcScenePrimitivesDescription("CexmcScenePrimitives")
void StackPreviousEvent(G4Event *anEvent)
void Start()
G4int ApplyCommand(const char *aCommand)
Definition: G4UImanager.cc:466
const G4String CexmcCalorimeterRegionName("Calorimeter")
void ProcessOneEvent(G4Event *anEvent)
#define G4endl
Definition: G4ios.hh:61
G4Timer * timer
static void PrintEnergyDeposit(const CexmcEnergyDepositStore *edStore)
static G4Eta * Definition()
Definition: G4Eta.cc:55
#define CEXMC_LINE_START
Definition: CexmcCommon.hh:52
G4UserEventAction * userEventAction
const G4String CexmcDetectorTypeName[CexmcNumberOfDetectorTypes]
G4bool EdTriggerIsOk(void) const
G4bool skipInteractionsWithoutEDTonWrite
G4Scene * GetCurrentScene() const
static G4UImanager * GetUIpointer()
Definition: G4UImanager.cc:73
G4HCofThisEvent * GetHCofThisEvent() const
Definition: G4Event.hh:185
static G4VVisManager * GetConcreteInstance()
G4double GetProposedMaxIL(void) const
CexmcEventDataVerboseLevel evDataVerboseLevel
virtual CexmcProductionModel * GetProductionModel(void)=0
G4bool IsFermiMotionOn(void) const
void SetProposedMaxIL(G4double value)
bool G4bool
Definition: G4Types.hh:79
CexmcProductionModelType productionModelType
static constexpr double deg
Definition: G4SIunits.hh:152
virtual G4Event * GenerateEvent(G4int i_event)
G4VUserPrimaryGeneratorAction * userPrimaryGeneratorAction
virtual void AnalyzeEvent(G4Event *anEvent)
void DumpInfo() const
G4int numberOfEventToBeProcessed
G4Run * currentRun
void Stop()
std::map< G4int, G4int > CexmcNmbOfHitsInRanges
Definition: CexmcRun.hh:51
G4bool ProjectIsSaved(void) const
G4bool ProjectIsRead(void) const
CexmcEventDataVerboseLevel rEvDataVerboseLevel
G4int numberOfEventsProcessedEffective
virtual ~CexmcRunManager()
G4bool runAborted
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
void SetAngularRanges(const CexmcAngularRangeList &angularRanges_)
CexmcRunManager(const G4String &projectId="", const G4String &rProject="", G4bool overrideExistingProject=false)
static G4SDManager * GetSDMpointer()
Definition: G4SDManager.cc:40
void ReadPreinitProjectData(void)
static G4RegionStore * GetInstance()
int G4int
Definition: G4Types.hh:78
std::vector< CexmcAngularRange > CexmcAngularRangeList
std::set< CexmcOutputDataType > CexmcOutputDataTypeSet
void RegisterScenePrimitives(void)
void SetupConstructionHook(void)
G4GLOB_DLL std::ostream G4cout
G4bool TpTriggerIsOk(void) const
void UpdateScoring()
static void PrintResults(const CexmcNmbOfHitsInRanges &nmbOfHitsSampled, const CexmcNmbOfHitsInRanges &nmbOfHitsSampledFull, const CexmcNmbOfHitsInRanges &nmbOfHitsTriggeredRealRange, const CexmcNmbOfHitsInRanges &nmbOfHitsTriggeredRecRange, const CexmcNmbOfHitsInRanges &nmbOfOrphanHits, const CexmcAngularRangeList &angularRanges, G4int nmbOfFalseHitsTriggeredEDT, G4int nmbOfFalseHitsTriggeredRec)
CexmcRunManagerMessenger * messenger
static void PrintTrackPoints(const CexmcTrackPointsStore *tpStore)
void BeamParticleChangeHook(void)
const CexmcAngularRangeList & GetAngularRanges(void) const
G4int verboseLevel
void ApplyFermiMotion(G4bool on, G4bool fromMessenger=true)
G4VUserEventInformation * GetUserInformation() const
Definition: G4Event.hh:199
CexmcEventCountPolicy eventCountPolicy
CexmcBasePhysicsUsed basePhysicsUsed
G4EventManager * eventManager