Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/examples/advanced/ChargeExchangeMC/src/CexmcEventAction.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /examples/advanced/ChargeExchangeMC/src/CexmcEventAction.cc (Version 11.3.0) and /examples/advanced/ChargeExchangeMC/src/CexmcEventAction.cc (Version 11.1)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                         3 // * License and Disclaimer                                           *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th      5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided      6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License      7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/      8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.           9 // * include a list of copyright holders.                             *
 10 // *                                               10 // *                                                                  *
 11 // * Neither the authors of this software syst     11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file      15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitatio     16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                               17 // *                                                                  *
 18 // * This  code  implementation is the result      18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboratio     19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distri     20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  ag     21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publicati     22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Sof     23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // *******************************************     24 // ********************************************************************
 25 //                                                 25 //
 26 /*                                                 26 /*
 27  * ===========================================     27  * ============================================================================
 28  *                                                 28  *
 29  *       Filename:  CexmcEventAction.cc            29  *       Filename:  CexmcEventAction.cc
 30  *                                                 30  *
 31  *    Description:  event action                   31  *    Description:  event action
 32  *                                                 32  *
 33  *        Version:  1.0                            33  *        Version:  1.0
 34  *        Created:  27.10.2009 22:48:08            34  *        Created:  27.10.2009 22:48:08
 35  *       Revision:  none                           35  *       Revision:  none
 36  *       Compiler:  gcc                            36  *       Compiler:  gcc
 37  *                                                 37  *
 38  *         Author:  Alexey Radkov (),              38  *         Author:  Alexey Radkov (), 
 39  *        Company:  PNPI                           39  *        Company:  PNPI
 40  *                                                 40  *
 41  * ===========================================     41  * ============================================================================
 42  */                                                42  */
 43                                                    43 
 44 #include <cmath>                                   44 #include <cmath>
 45 #ifdef CEXMC_USE_PERSISTENCY                       45 #ifdef CEXMC_USE_PERSISTENCY
 46 #include <boost/archive/binary_oarchive.hpp>       46 #include <boost/archive/binary_oarchive.hpp>
 47 #endif                                             47 #endif
 48 #include <G4DigiManager.hh>                        48 #include <G4DigiManager.hh>
 49 #include <G4Event.hh>                              49 #include <G4Event.hh>
 50 #include <G4Circle.hh>                             50 #include <G4Circle.hh>
 51 #include <G4VisAttributes.hh>                      51 #include <G4VisAttributes.hh>
 52 #include <G4VisManager.hh>                         52 #include <G4VisManager.hh>
 53 #include <G4VTrajectory.hh>                        53 #include <G4VTrajectory.hh>
 54 #include <G4TrajectoryContainer.hh>                54 #include <G4TrajectoryContainer.hh>
 55 #include <G4Colour.hh>                             55 #include <G4Colour.hh>
 56 #include <G4SystemOfUnits.hh>                      56 #include <G4SystemOfUnits.hh>
 57 #include "CexmcEventAction.hh"                     57 #include "CexmcEventAction.hh"
 58 #include "CexmcEventActionMessenger.hh"            58 #include "CexmcEventActionMessenger.hh"
 59 #include "CexmcEventInfo.hh"                       59 #include "CexmcEventInfo.hh"
 60 #include "CexmcEventSObject.hh"                    60 #include "CexmcEventSObject.hh"
 61 #include "CexmcEventFastSObject.hh"                61 #include "CexmcEventFastSObject.hh"
 62 #include "CexmcTrackingAction.hh"                  62 #include "CexmcTrackingAction.hh"
 63 #include "CexmcChargeExchangeReconstructor.hh"     63 #include "CexmcChargeExchangeReconstructor.hh"
 64 #include "CexmcRunManager.hh"                      64 #include "CexmcRunManager.hh"
 65 #include "CexmcHistoManager.hh"                    65 #include "CexmcHistoManager.hh"
 66 #include "CexmcRun.hh"                             66 #include "CexmcRun.hh"
 67 #include "CexmcPhysicsManager.hh"                  67 #include "CexmcPhysicsManager.hh"
 68 #include "CexmcProductionModel.hh"                 68 #include "CexmcProductionModel.hh"
 69 #include "CexmcProductionModelData.hh"             69 #include "CexmcProductionModelData.hh"
 70 #include "CexmcEnergyDepositDigitizer.hh"          70 #include "CexmcEnergyDepositDigitizer.hh"
 71 #include "CexmcEnergyDepositStore.hh"              71 #include "CexmcEnergyDepositStore.hh"
 72 #include "CexmcTrackPointsDigitizer.hh"            72 #include "CexmcTrackPointsDigitizer.hh"
 73 #include "CexmcTrackPointsStore.hh"                73 #include "CexmcTrackPointsStore.hh"
 74 #include "CexmcTrackPointInfo.hh"                  74 #include "CexmcTrackPointInfo.hh"
 75 #include "CexmcException.hh"                       75 #include "CexmcException.hh"
 76 #include "CexmcCommon.hh"                          76 #include "CexmcCommon.hh"
 77                                                    77 
 78                                                    78 
 79 namespace                                          79 namespace
 80 {                                                  80 {
 81     G4double  CexmcSmallCircleScreenSize( 5.0      81     G4double  CexmcSmallCircleScreenSize( 5.0 );
 82     G4double  CexmcBigCircleScreenSize( 10.0 )     82     G4double  CexmcBigCircleScreenSize( 10.0 );
 83     G4Colour  CexmcTrackPointsMarkerColour( 0.     83     G4Colour  CexmcTrackPointsMarkerColour( 0.0, 1.0, 0.4 );
 84     G4Colour  CexmcRecTrackPointsMarkerColour(     84     G4Colour  CexmcRecTrackPointsMarkerColour( 1.0, 0.4, 0.0 );
 85                                                    85 
 86 #ifdef CEXMC_USE_ROOT                              86 #ifdef CEXMC_USE_ROOT
 87     inline G4double  CexmcGetKinEnergy( G4doub     87     inline G4double  CexmcGetKinEnergy( G4double  momentumAmp, G4double  mass )
 88     {                                              88     {
 89         return std::sqrt( momentumAmp * moment     89         return std::sqrt( momentumAmp * momentumAmp + mass * mass ) - mass;
 90     }                                              90     }
 91 #endif                                             91 #endif
 92 }                                                  92 }
 93                                                    93 
 94                                                    94 
 95 CexmcEventAction::CexmcEventAction( CexmcPhysi     95 CexmcEventAction::CexmcEventAction( CexmcPhysicsManager *  physicsManager_,
 96                                     G4int  ver     96                                     G4int  verbose_ ) :
 97     physicsManager( physicsManager_ ), reconst     97     physicsManager( physicsManager_ ), reconstructor( NULL ),
 98 #ifdef CEXMC_USE_ROOT                              98 #ifdef CEXMC_USE_ROOT
 99     opKinEnergy( 0. ),                             99     opKinEnergy( 0. ),
100 #endif                                            100 #endif
101     verbose( verbose_ ), verboseDraw( 4 ), mes    101     verbose( verbose_ ), verboseDraw( 4 ), messenger( NULL )
102 {                                                 102 {
103     G4DigiManager *  digiManager( G4DigiManage    103     G4DigiManager *  digiManager( G4DigiManager::GetDMpointer() );
104     digiManager->AddNewModule( new CexmcEnergy    104     digiManager->AddNewModule( new CexmcEnergyDepositDigitizer(
105                                                   105                                                     CexmcEDDigitizerName ) );
106     digiManager->AddNewModule( new CexmcTrackP    106     digiManager->AddNewModule( new CexmcTrackPointsDigitizer(
107                                                   107                                                     CexmcTPDigitizerName ) );
108     reconstructor = new CexmcChargeExchangeRec    108     reconstructor = new CexmcChargeExchangeReconstructor(
109                                         physic    109                                         physicsManager->GetProductionModel() );
110     messenger = new CexmcEventActionMessenger(    110     messenger = new CexmcEventActionMessenger( this );
111 }                                                 111 }
112                                                   112 
113                                                   113 
114 CexmcEventAction::~CexmcEventAction()             114 CexmcEventAction::~CexmcEventAction()
115 {                                                 115 {
116     delete reconstructor;                         116     delete reconstructor;
117     delete messenger;                             117     delete messenger;
118 }                                                 118 }
119                                                   119 
120                                                   120 
121 void  CexmcEventAction::BeamParticleChangeHook    121 void  CexmcEventAction::BeamParticleChangeHook( void )
122 {                                                 122 {
123     reconstructor->SetupBeamParticle();           123     reconstructor->SetupBeamParticle();
124 }                                                 124 }
125                                                   125 
126                                                   126 
127 void  CexmcEventAction::BeginOfEventAction( co    127 void  CexmcEventAction::BeginOfEventAction( const G4Event * )
128 {                                                 128 {
129     G4RunManager *         runManager( G4RunMa    129     G4RunManager *         runManager( G4RunManager::GetRunManager() );
130     CexmcTrackingAction *  trackingAction         130     CexmcTrackingAction *  trackingAction
131             ( static_cast< CexmcTrackingAction    131             ( static_cast< CexmcTrackingAction * >(
132                         const_cast< G4UserTrac    132                         const_cast< G4UserTrackingAction * >(
133                                     runManager    133                                     runManager->GetUserTrackingAction() ) ) );
134     trackingAction->BeginOfEventAction();         134     trackingAction->BeginOfEventAction();
135                                                   135 
136     physicsManager->ResetNumberOfTriggeredStud    136     physicsManager->ResetNumberOfTriggeredStudiedInteractions();
137 }                                                 137 }
138                                                   138 
139                                                   139 
140 CexmcEnergyDepositStore *  CexmcEventAction::M    140 CexmcEnergyDepositStore *  CexmcEventAction::MakeEnergyDepositStore(
141                                 const CexmcEne    141                                 const CexmcEnergyDepositDigitizer *  digitizer )
142 {                                                 142 {
143     G4double  monitorED( digitizer->GetMonitor    143     G4double  monitorED( digitizer->GetMonitorED() );
144     G4double  vetoCounterEDLeft( digitizer->Ge    144     G4double  vetoCounterEDLeft( digitizer->GetVetoCounterEDLeft() );
145     G4double  vetoCounterEDRight( digitizer->G    145     G4double  vetoCounterEDRight( digitizer->GetVetoCounterEDRight() );
146     G4double  calorimeterEDLeft( digitizer->Ge    146     G4double  calorimeterEDLeft( digitizer->GetCalorimeterEDLeft() );
147     G4double  calorimeterEDRight( digitizer->G    147     G4double  calorimeterEDRight( digitizer->GetCalorimeterEDRight() );
148     G4int     calorimeterEDLeftMaxX( digitizer    148     G4int     calorimeterEDLeftMaxX( digitizer->GetCalorimeterEDLeftMaxX() );
149     G4int     calorimeterEDLeftMaxY( digitizer    149     G4int     calorimeterEDLeftMaxY( digitizer->GetCalorimeterEDLeftMaxY() );
150     G4int     calorimeterEDRightMaxX( digitize    150     G4int     calorimeterEDRightMaxX( digitizer->GetCalorimeterEDRightMaxX() );
151     G4int     calorimeterEDRightMaxY( digitize    151     G4int     calorimeterEDRightMaxY( digitizer->GetCalorimeterEDRightMaxY() );
152                                                   152 
153     const CexmcEnergyDepositCalorimeterCollect    153     const CexmcEnergyDepositCalorimeterCollection &
154               calorimeterEDLeftCollection(        154               calorimeterEDLeftCollection(
155                             digitizer->GetCalo    155                             digitizer->GetCalorimeterEDLeftCollection() );
156     const CexmcEnergyDepositCalorimeterCollect    156     const CexmcEnergyDepositCalorimeterCollection &
157               calorimeterEDRightCollection(       157               calorimeterEDRightCollection(
158                             digitizer->GetCalo    158                             digitizer->GetCalorimeterEDRightCollection() );
159                                                   159 
160     /* ATTENTION: return object in heap - must    160     /* ATTENTION: return object in heap - must be freed by caller! */
161     return new CexmcEnergyDepositStore( monito    161     return new CexmcEnergyDepositStore( monitorED, vetoCounterEDLeft,
162                     vetoCounterEDRight, calori    162                     vetoCounterEDRight, calorimeterEDLeft, calorimeterEDRight,
163                     calorimeterEDLeftMaxX, cal    163                     calorimeterEDLeftMaxX, calorimeterEDLeftMaxY,
164                     calorimeterEDRightMaxX, ca    164                     calorimeterEDRightMaxX, calorimeterEDRightMaxY,
165                     calorimeterEDLeftCollectio    165                     calorimeterEDLeftCollection, calorimeterEDRightCollection );
166 }                                                 166 }
167                                                   167 
168                                                   168 
169 CexmcTrackPointsStore *  CexmcEventAction::Mak    169 CexmcTrackPointsStore *  CexmcEventAction::MakeTrackPointsStore(
170                                 const CexmcTra    170                                 const CexmcTrackPointsDigitizer *  digitizer )
171 {                                                 171 {
172     const CexmcTrackPointInfo &                   172     const CexmcTrackPointInfo &
173                 monitorTP( digitizer->GetMonit    173                 monitorTP( digitizer->GetMonitorTP() );
174     const CexmcTrackPointInfo &                   174     const CexmcTrackPointInfo &
175                 targetTPBeamParticle(             175                 targetTPBeamParticle(
176                     digitizer->GetTargetTPBeam    176                     digitizer->GetTargetTPBeamParticle() );
177     const CexmcTrackPointInfo &                   177     const CexmcTrackPointInfo &
178                 targetTPOutputParticle(           178                 targetTPOutputParticle(
179                     digitizer->GetTargetTPOutp    179                     digitizer->GetTargetTPOutputParticle() );
180     const CexmcTrackPointInfo &                   180     const CexmcTrackPointInfo &
181                 targetTPNucleusParticle(          181                 targetTPNucleusParticle(
182                     digitizer->GetTargetTPNucl    182                     digitizer->GetTargetTPNucleusParticle() );
183     const CexmcTrackPointInfo &                   183     const CexmcTrackPointInfo &
184                 targetTPOutputParticleDecayPro    184                 targetTPOutputParticleDecayProductParticle1(
185                     digitizer->                   185                     digitizer->
186                     GetTargetTPOutputParticleD    186                     GetTargetTPOutputParticleDecayProductParticle( 0 ) );
187     const CexmcTrackPointInfo &                   187     const CexmcTrackPointInfo &
188                 targetTPOutputParticleDecayPro    188                 targetTPOutputParticleDecayProductParticle2(
189                     digitizer->                   189                     digitizer->
190                     GetTargetTPOutputParticleD    190                     GetTargetTPOutputParticleDecayProductParticle( 1 ) );
191     const CexmcTrackPointInfo &                   191     const CexmcTrackPointInfo &
192                 vetoCounterTPLeft(                192                 vetoCounterTPLeft(
193                     digitizer->GetVetoCounterT    193                     digitizer->GetVetoCounterTPLeft() );
194     const CexmcTrackPointInfo &                   194     const CexmcTrackPointInfo &
195                 vetoCounterTPRight(               195                 vetoCounterTPRight(
196                     digitizer->GetVetoCounterT    196                     digitizer->GetVetoCounterTPRight() );
197     const CexmcTrackPointInfo &                   197     const CexmcTrackPointInfo &
198                 calorimeterTPLeft(                198                 calorimeterTPLeft(
199                     digitizer->GetCalorimeterT    199                     digitizer->GetCalorimeterTPLeft() );
200     const CexmcTrackPointInfo &                   200     const CexmcTrackPointInfo &
201                 calorimeterTPRight(               201                 calorimeterTPRight(
202                     digitizer->GetCalorimeterT    202                     digitizer->GetCalorimeterTPRight() );
203                                                   203 
204     /* ATTENTION: return object in heap - must    204     /* ATTENTION: return object in heap - must be freed by caller! */
205     return new CexmcTrackPointsStore( monitorT    205     return new CexmcTrackPointsStore( monitorTP, targetTPBeamParticle,
206                   targetTPOutputParticle, targ    206                   targetTPOutputParticle, targetTPNucleusParticle,
207                   targetTPOutputParticleDecayP    207                   targetTPOutputParticleDecayProductParticle1,
208                   targetTPOutputParticleDecayP    208                   targetTPOutputParticleDecayProductParticle2,
209                   vetoCounterTPLeft, vetoCount    209                   vetoCounterTPLeft, vetoCounterTPRight,
210                   calorimeterTPLeft, calorimet    210                   calorimeterTPLeft, calorimeterTPRight );
211 }                                                 211 }
212                                                   212 
213                                                   213 
214 void  CexmcEventAction::PrintEnergyDeposit(       214 void  CexmcEventAction::PrintEnergyDeposit(
215                                     const Cexm    215                                     const CexmcEnergyDepositStore *  edStore )
216 {                                                 216 {
217     G4cout << " --- Energy Deposit" << G4endl;    217     G4cout << " --- Energy Deposit" << G4endl;
218     G4cout << "       monitor : " <<              218     G4cout << "       monitor : " <<
219             G4BestUnit( edStore->monitorED, "E    219             G4BestUnit( edStore->monitorED, "Energy" ) << G4endl;
220     G4cout << "        vc (l) : " <<              220     G4cout << "        vc (l) : " <<
221             G4BestUnit( edStore->vetoCounterED    221             G4BestUnit( edStore->vetoCounterEDLeft, "Energy" ) << G4endl;
222     G4cout << "        vc (r) : " <<              222     G4cout << "        vc (r) : " <<
223             G4BestUnit( edStore->vetoCounterED    223             G4BestUnit( edStore->vetoCounterEDRight, "Energy" ) << G4endl;
224     G4cout << "       cal (l) : " <<              224     G4cout << "       cal (l) : " <<
225             G4BestUnit( edStore->calorimeterED    225             G4BestUnit( edStore->calorimeterEDLeft, "Energy" );
226     G4cout << edStore->calorimeterEDLeftCollec    226     G4cout << edStore->calorimeterEDLeftCollection;
227     G4cout << "       cal (r) : " <<              227     G4cout << "       cal (r) : " <<
228             G4BestUnit( edStore->calorimeterED    228             G4BestUnit( edStore->calorimeterEDRight, "Energy" );
229     G4cout << edStore->calorimeterEDRightColle    229     G4cout << edStore->calorimeterEDRightCollection;
230 }                                                 230 }
231                                                   231 
232                                                   232 
233 void  CexmcEventAction::PrintTrackPoints(         233 void  CexmcEventAction::PrintTrackPoints(
234                                     const Cexm    234                                     const CexmcTrackPointsStore *  tpStore )
235 {                                                 235 {
236     if ( ! tpStore )                              236     if ( ! tpStore )
237         return;                                   237         return;
238                                                   238 
239     G4cout << " --- Track Points" << G4endl;      239     G4cout << " --- Track Points" << G4endl;
240     G4cout << "       monitor : " << tpStore->    240     G4cout << "       monitor : " << tpStore->monitorTP << G4endl;
241     G4cout << "        target : " << tpStore->    241     G4cout << "        target : " << tpStore->targetTPBeamParticle << G4endl;
242     G4cout << "               : " << tpStore->    242     G4cout << "               : " << tpStore->targetTPOutputParticle << G4endl;
243     G4cout << "               : " << tpStore->    243     G4cout << "               : " << tpStore->targetTPNucleusParticle << G4endl;
244     G4cout << "               : " <<              244     G4cout << "               : " <<
245             tpStore->targetTPOutputParticleDec    245             tpStore->targetTPOutputParticleDecayProductParticle1 << G4endl;
246     G4cout << "               : " <<              246     G4cout << "               : " <<
247             tpStore->targetTPOutputParticleDec    247             tpStore->targetTPOutputParticleDecayProductParticle2 << G4endl;
248     G4cout << "        vc (l) : " << tpStore->    248     G4cout << "        vc (l) : " << tpStore->vetoCounterTPLeft << G4endl;
249     G4cout << "        vc (r) : " << tpStore->    249     G4cout << "        vc (r) : " << tpStore->vetoCounterTPRight << G4endl;
250     G4cout << "       cal (l) : " << tpStore->    250     G4cout << "       cal (l) : " << tpStore->calorimeterTPLeft << G4endl;
251     G4cout << "       cal (r) : " << tpStore->    251     G4cout << "       cal (r) : " << tpStore->calorimeterTPRight << G4endl;
252     G4cout << "      ---" << G4endl;              252     G4cout << "      ---" << G4endl;
253     G4cout << "      angle between the " <<       253     G4cout << "      angle between the " <<
254         tpStore->targetTPOutputParticle.partic    254         tpStore->targetTPOutputParticle.particle->GetParticleName() <<
255         " decay products : " <<                   255         " decay products : " <<
256         tpStore->targetTPOutputParticleDecayPr    256         tpStore->targetTPOutputParticleDecayProductParticle1.directionWorld.
257             angle(                                257             angle(
258         tpStore->targetTPOutputParticleDecayPr    258         tpStore->targetTPOutputParticleDecayProductParticle2.directionWorld ) /
259             deg << " deg" << G4endl;              259             deg << " deg" << G4endl;
260 }                                                 260 }
261                                                   261 
262                                                   262 
263 void  CexmcEventAction::PrintProductionModelDa    263 void  CexmcEventAction::PrintProductionModelData(
264                                 const CexmcAng    264                                 const CexmcAngularRangeList &  angularRanges,
265                                 const CexmcPro    265                                 const CexmcProductionModelData &  pmData )
266 {                                                 266 {
267     G4cout << " --- Triggered angular ranges:     267     G4cout << " --- Triggered angular ranges: " << angularRanges;
268     G4cout << " --- Production model data: " <    268     G4cout << " --- Production model data: " << pmData;
269 }                                                 269 }
270                                                   270 
271                                                   271 
272 void  CexmcEventAction::PrintReconstructedData    272 void  CexmcEventAction::PrintReconstructedData(
273                     const CexmcAngularRangeLis    273                     const CexmcAngularRangeList &  triggeredRecAngularRanges,
274                     const CexmcAngularRange &     274                     const CexmcAngularRange &  angularGap ) const
275 {                                                 275 {
276     G4cout << " --- Reconstructed data: " << G    276     G4cout << " --- Reconstructed data: " << G4endl;
277     G4cout << "       -- entry points:" << G4e    277     G4cout << "       -- entry points:" << G4endl;
278     G4cout << "              left: " << G4Best    278     G4cout << "              left: " << G4BestUnit(
279         reconstructor->GetCalorimeterEPLeftPos    279         reconstructor->GetCalorimeterEPLeftPosition(), "Length" ) << G4endl;
280     G4cout << "             right: " << G4Best    280     G4cout << "             right: " << G4BestUnit(
281         reconstructor->GetCalorimeterEPRightPo    281         reconstructor->GetCalorimeterEPRightPosition(), "Length" ) << G4endl;
282     G4cout << "            target: " << G4Best    282     G4cout << "            target: " << G4BestUnit(
283         reconstructor->GetTargetEPPosition(),     283         reconstructor->GetTargetEPPosition(), "Length" ) << G4endl;
284     G4cout << "       -- the angle: " << recon    284     G4cout << "       -- the angle: " << reconstructor->GetTheAngle() / deg <<
285         " deg" << G4endl;                         285         " deg" << G4endl;
286     G4cout << "       -- mass of the output pa    286     G4cout << "       -- mass of the output particle: " << G4BestUnit(
287         reconstructor->GetOutputParticleMass()    287         reconstructor->GetOutputParticleMass(), "Energy" ) << G4endl;
288     G4cout << "       -- mass of the nucleus o    288     G4cout << "       -- mass of the nucleus output particle: " << G4BestUnit(
289         reconstructor->GetNucleusOutputParticl    289         reconstructor->GetNucleusOutputParticleMass(), "Energy" ) << G4endl;
290     if ( reconstructor->IsMassCutUsed() )         290     if ( reconstructor->IsMassCutUsed() )
291     {                                             291     {
292         if ( reconstructor->HasMassCutTriggere    292         if ( reconstructor->HasMassCutTriggered() )
293             G4cout << "            < mass cut     293             G4cout << "            < mass cut passed >" << G4endl;
294         else                                      294         else
295             G4cout << "            < mass cut     295             G4cout << "            < mass cut failed >" << G4endl;
296     }                                             296     }
297     if ( reconstructor->IsAbsorbedEnergyCutUse    297     if ( reconstructor->IsAbsorbedEnergyCutUsed() )
298     {                                             298     {
299         if ( reconstructor->HasAbsorbedEnergyC    299         if ( reconstructor->HasAbsorbedEnergyCutTriggered() )
300             G4cout << "            < absorbed     300             G4cout << "            < absorbed energy cut passed >" << G4endl;
301         else                                      301         else
302             G4cout << "            < absorbed     302             G4cout << "            < absorbed energy cut failed >" << G4endl;
303     }                                             303     }
304     const CexmcProductionModelData &  pmData(     304     const CexmcProductionModelData &  pmData(
305                                     reconstruc    305                                     reconstructor->GetProductionModelData() );
306     G4cout << "       -- production model data    306     G4cout << "       -- production model data: " << pmData;
307     G4cout << "       -- triggered angular ran    307     G4cout << "       -- triggered angular ranges: ";
308     if ( triggeredRecAngularRanges.empty() )      308     if ( triggeredRecAngularRanges.empty() )
309         G4cout << "< orphan detected, gap " <<    309         G4cout << "< orphan detected, gap " << angularGap << " >" << G4endl;
310     else                                          310     else
311         G4cout << triggeredRecAngularRanges;      311         G4cout << triggeredRecAngularRanges;
312 }                                                 312 }
313                                                   313 
314                                                   314 
315 #ifdef CEXMC_USE_ROOT                             315 #ifdef CEXMC_USE_ROOT
316                                                   316 
317 void  CexmcEventAction::FillEDTHistos( const C    317 void  CexmcEventAction::FillEDTHistos( const CexmcEnergyDepositStore *  edStore,
318                 const CexmcAngularRangeList &     318                 const CexmcAngularRangeList &  triggeredAngularRanges ) const
319 {                                                 319 {
320     CexmcHistoManager *  histoManager( CexmcHi    320     CexmcHistoManager *  histoManager( CexmcHistoManager::Instance() );
321                                                   321 
322     histoManager->Add( CexmcAbsorbedEnergy_EDT    322     histoManager->Add( CexmcAbsorbedEnergy_EDT_Histo, 0,
323                        edStore->calorimeterEDL    323                        edStore->calorimeterEDLeft,
324                        edStore->calorimeterEDR    324                        edStore->calorimeterEDRight );
325                                                   325 
326     for ( CexmcAngularRangeList::const_iterato    326     for ( CexmcAngularRangeList::const_iterator
327                         k( triggeredAngularRan    327                         k( triggeredAngularRanges.begin() );
328                                     k != trigg    328                                     k != triggeredAngularRanges.end(); ++k )
329     {                                             329     {
330         histoManager->Add( CexmcAbsEnInLeftCal    330         histoManager->Add( CexmcAbsEnInLeftCalorimeter_ARReal_EDT_Histo,
331                            k->index, edStore->    331                            k->index, edStore->calorimeterEDLeft );
332         histoManager->Add( CexmcAbsEnInRightCa    332         histoManager->Add( CexmcAbsEnInRightCalorimeter_ARReal_EDT_Histo,
333                            k->index, edStore->    333                            k->index, edStore->calorimeterEDRight );
334     }                                             334     }
335 }                                                 335 }
336                                                   336 
337                                                   337 
338 void  CexmcEventAction::FillTPTHistos( const C    338 void  CexmcEventAction::FillTPTHistos( const CexmcTrackPointsStore *  tpStore,
339                 const CexmcProductionModelData    339                 const CexmcProductionModelData &  pmData,
340                 const CexmcAngularRangeList &     340                 const CexmcAngularRangeList &  triggeredAngularRanges ) const
341 {                                                 341 {
342     CexmcHistoManager *  histoManager( CexmcHi    342     CexmcHistoManager *  histoManager( CexmcHistoManager::Instance() );
343                                                   343 
344     if ( tpStore->monitorTP.IsValid() )           344     if ( tpStore->monitorTP.IsValid() )
345     {                                             345     {
346         histoManager->Add( CexmcMomentumBP_TPT    346         histoManager->Add( CexmcMomentumBP_TPT_Histo, 0,
347                            tpStore->monitorTP.    347                            tpStore->monitorTP.momentumAmp );
348         histoManager->Add( CexmcTPInMonitor_TP    348         histoManager->Add( CexmcTPInMonitor_TPT_Histo, 0,
349                            tpStore->monitorTP.    349                            tpStore->monitorTP.positionLocal.x(),
350                            tpStore->monitorTP.    350                            tpStore->monitorTP.positionLocal.y() );
351     }                                             351     }
352                                                   352 
353     if ( tpStore->targetTPOutputParticle.IsVal    353     if ( tpStore->targetTPOutputParticle.IsValid() )
354     {                                             354     {
355         histoManager->Add( CexmcTPInTarget_TPT    355         histoManager->Add( CexmcTPInTarget_TPT_Histo, 0,
356                            tpStore->targetTPOu    356                            tpStore->targetTPOutputParticle.positionLocal.x(),
357                            tpStore->targetTPOu    357                            tpStore->targetTPOutputParticle.positionLocal.y(),
358                            tpStore->targetTPOu    358                            tpStore->targetTPOutputParticle.positionLocal.z() );
359         if ( histoManager->GetVerboseLevel() >    359         if ( histoManager->GetVerboseLevel() > 0 )
360         {                                         360         {
361             histoManager->Add( CexmcMomentumIP    361             histoManager->Add( CexmcMomentumIP_TPT_Histo, 0,
362                                pmData.incident    362                                pmData.incidentParticleLAB.rho() );
363         }                                         363         }
364     }                                             364     }
365                                                   365 
366     for ( CexmcAngularRangeList::const_iterato    366     for ( CexmcAngularRangeList::const_iterator
367                         k( triggeredAngularRan    367                         k( triggeredAngularRanges.begin() );
368                                     k != trigg    368                                     k != triggeredAngularRanges.end(); ++k )
369     {                                             369     {
370         if ( tpStore->calorimeterTPLeft.IsVali    370         if ( tpStore->calorimeterTPLeft.IsValid() )
371         {                                         371         {
372             /* kinetic energy and momentum of     372             /* kinetic energy and momentum of gamma are equal */
373             G4double  kinEnergy( tpStore->calo    373             G4double  kinEnergy( tpStore->calorimeterTPLeft.momentumAmp );
374             histoManager->Add( CexmcKinEnAtLef    374             histoManager->Add( CexmcKinEnAtLeftCalorimeter_ARReal_TPT_Histo,
375                                k->index, kinEn    375                                k->index, kinEnergy );
376         }                                         376         }
377         if ( tpStore->calorimeterTPRight.IsVal    377         if ( tpStore->calorimeterTPRight.IsValid() )
378         {                                         378         {
379             G4double  kinEnergy( tpStore->calo    379             G4double  kinEnergy( tpStore->calorimeterTPRight.momentumAmp );
380             histoManager->Add( CexmcKinEnAtRig    380             histoManager->Add( CexmcKinEnAtRightCalorimeter_ARReal_TPT_Histo,
381                                k->index, kinEn    381                                k->index, kinEnergy );
382         }                                         382         }
383         if ( tpStore->targetTPOutputParticle.I    383         if ( tpStore->targetTPOutputParticle.IsValid() )
384         {                                         384         {
385             histoManager->Add( CexmcTPInTarget    385             histoManager->Add( CexmcTPInTarget_ARReal_TPT_Histo, k->index,
386                            tpStore->targetTPOu    386                            tpStore->targetTPOutputParticle.positionLocal.x(),
387                            tpStore->targetTPOu    387                            tpStore->targetTPOutputParticle.positionLocal.y(),
388                            tpStore->targetTPOu    388                            tpStore->targetTPOutputParticle.positionLocal.z() );
389             histoManager->Add( CexmcKinEnOP_LA    389             histoManager->Add( CexmcKinEnOP_LAB_ARReal_TPT_Histo, k->index,
390                                opKinEnergy );     390                                opKinEnergy );
391             histoManager->Add( CexmcAngleOP_SC    391             histoManager->Add( CexmcAngleOP_SCM_ARReal_TPT_Histo, k->index,
392                                pmData.outputPa    392                                pmData.outputParticleSCM.cosTheta() );
393         }                                         393         }
394         if ( tpStore->targetTPOutputParticleDe    394         if ( tpStore->targetTPOutputParticleDecayProductParticle1.IsValid() &&
395              tpStore->targetTPOutputParticleDe    395              tpStore->targetTPOutputParticleDecayProductParticle2.IsValid() )
396         {                                         396         {
397             G4double  openAngle(                  397             G4double  openAngle(
398                         tpStore->targetTPOutpu    398                         tpStore->targetTPOutputParticleDecayProductParticle1.
399                         directionWorld.angle(     399                         directionWorld.angle( tpStore->
400                                 targetTPOutput    400                                 targetTPOutputParticleDecayProductParticle2.
401                                         direct    401                                         directionWorld ) / deg );
402             histoManager->Add( CexmcOpenAngle_    402             histoManager->Add( CexmcOpenAngle_ARReal_TPT_Histo, k->index,
403                                openAngle );       403                                openAngle );
404         }                                         404         }
405     }                                             405     }
406 }                                                 406 }
407                                                   407 
408                                                   408 
409 void  CexmcEventAction::FillRTHistos( G4bool      409 void  CexmcEventAction::FillRTHistos( G4bool  reconstructorHasFullTrigger,
410                 const CexmcEnergyDepositStore     410                 const CexmcEnergyDepositStore *  edStore,
411                 const CexmcTrackPointsStore *     411                 const CexmcTrackPointsStore *  tpStore,
412                 const CexmcProductionModelData    412                 const CexmcProductionModelData &  pmData,
413                 const CexmcAngularRangeList &     413                 const CexmcAngularRangeList &  triggeredAngularRanges ) const
414 {                                                 414 {
415     CexmcHistoManager *  histoManager( CexmcHi    415     CexmcHistoManager *  histoManager( CexmcHistoManager::Instance() );
416                                                   416 
417     G4double    opMass( reconstructor->GetOutp    417     G4double    opMass( reconstructor->GetOutputParticleMass() );
418     G4double    nopMass( reconstructor->GetNuc    418     G4double    nopMass( reconstructor->GetNucleusOutputParticleMass() );
419                                                   419 
420     histoManager->Add( CexmcRecMasses_EDT_Hist    420     histoManager->Add( CexmcRecMasses_EDT_Histo, 0, opMass, nopMass );
421                                                   421 
422     for ( CexmcAngularRangeList::const_iterato    422     for ( CexmcAngularRangeList::const_iterator
423                         k( triggeredAngularRan    423                         k( triggeredAngularRanges.begin() );
424                                     k != trigg    424                                     k != triggeredAngularRanges.end(); ++k )
425     {                                             425     {
426         if ( tpStore->calorimeterTPLeft.IsVali    426         if ( tpStore->calorimeterTPLeft.IsValid() )
427         {                                         427         {
428             histoManager->Add( CexmcOPDPAtLeft    428             histoManager->Add( CexmcOPDPAtLeftCalorimeter_ARReal_EDT_Histo,
429                                k->index,          429                                k->index,
430                                tpStore->calori    430                                tpStore->calorimeterTPLeft.positionLocal.x(),
431                                tpStore->calori    431                                tpStore->calorimeterTPLeft.positionLocal.y() );
432         }                                         432         }
433         if ( tpStore->calorimeterTPRight.IsVal    433         if ( tpStore->calorimeterTPRight.IsValid() )
434         {                                         434         {
435             histoManager->Add( CexmcOPDPAtRigh    435             histoManager->Add( CexmcOPDPAtRightCalorimeter_ARReal_EDT_Histo,
436                                k->index,          436                                k->index,
437                                tpStore->calori    437                                tpStore->calorimeterTPRight.positionLocal.x(),
438                                tpStore->calori    438                                tpStore->calorimeterTPRight.positionLocal.y() );
439         }                                         439         }
440         histoManager->Add( CexmcRecOPDPAtLeftC    440         histoManager->Add( CexmcRecOPDPAtLeftCalorimeter_ARReal_EDT_Histo,
441                            k->index,              441                            k->index,
442                            reconstructor->GetC    442                            reconstructor->GetCalorimeterEPLeftPosition().x(),
443                            reconstructor->GetC    443                            reconstructor->GetCalorimeterEPLeftPosition().y() );
444         histoManager->Add( CexmcRecOPDPAtRight    444         histoManager->Add( CexmcRecOPDPAtRightCalorimeter_ARReal_EDT_Histo,
445                            k->index,              445                            k->index,
446                            reconstructor->GetC    446                            reconstructor->GetCalorimeterEPRightPosition().x(),
447                            reconstructor->GetC    447                            reconstructor->GetCalorimeterEPRightPosition().y() );
448     }                                             448     }
449                                                   449 
450     if ( ! reconstructorHasFullTrigger )          450     if ( ! reconstructorHasFullTrigger )
451         return;                                   451         return;
452                                                   452 
453     if ( tpStore->monitorTP.IsValid() )           453     if ( tpStore->monitorTP.IsValid() )
454     {                                             454     {
455         histoManager->Add( CexmcMomentumBP_RT_    455         histoManager->Add( CexmcMomentumBP_RT_Histo, 0,
456                            tpStore->monitorTP.    456                            tpStore->monitorTP.momentumAmp );
457     }                                             457     }
458                                                   458 
459     if ( tpStore->targetTPOutputParticle.IsVal    459     if ( tpStore->targetTPOutputParticle.IsValid() )
460     {                                             460     {
461         histoManager->Add( CexmcTPInTarget_RT_    461         histoManager->Add( CexmcTPInTarget_RT_Histo, 0,
462                            tpStore->targetTPOu    462                            tpStore->targetTPOutputParticle.positionLocal.x(),
463                            tpStore->targetTPOu    463                            tpStore->targetTPOutputParticle.positionLocal.y(),
464                            tpStore->targetTPOu    464                            tpStore->targetTPOutputParticle.positionLocal.z() );
465     }                                             465     }
466                                                   466 
467     histoManager->Add( CexmcRecMasses_RT_Histo    467     histoManager->Add( CexmcRecMasses_RT_Histo, 0,
468                        reconstructor->GetOutpu    468                        reconstructor->GetOutputParticleMass(),
469                        reconstructor->GetNucle    469                        reconstructor->GetNucleusOutputParticleMass() );
470                                                   470 
471     histoManager->Add( CexmcAbsorbedEnergy_RT_    471     histoManager->Add( CexmcAbsorbedEnergy_RT_Histo, 0,
472                        edStore->calorimeterEDL    472                        edStore->calorimeterEDLeft,
473                        edStore->calorimeterEDR    473                        edStore->calorimeterEDRight );
474                                                   474 
475     G4double  recCosTheta( reconstructor->GetP    475     G4double  recCosTheta( reconstructor->GetProductionModelData().
476                                                   476                                                outputParticleSCM.cosTheta());
477                                                   477 
478     for ( CexmcAngularRangeList::const_iterato    478     for ( CexmcAngularRangeList::const_iterator
479                         k( triggeredAngularRan    479                         k( triggeredAngularRanges.begin() );
480                                     k != trigg    480                                     k != triggeredAngularRanges.end(); ++k )
481     {                                             481     {
482         histoManager->Add( CexmcRecMassOP_ARRe    482         histoManager->Add( CexmcRecMassOP_ARReal_RT_Histo, k->index, opMass );
483         histoManager->Add( CexmcRecMassNOP_ARR    483         histoManager->Add( CexmcRecMassNOP_ARReal_RT_Histo, k->index, nopMass );
484         if ( tpStore->calorimeterTPLeft.IsVali    484         if ( tpStore->calorimeterTPLeft.IsValid() )
485         {                                         485         {
486             G4double  kinEnergy( tpStore->calo    486             G4double  kinEnergy( tpStore->calorimeterTPLeft.momentumAmp );
487             histoManager->Add( CexmcKinEnAtLef    487             histoManager->Add( CexmcKinEnAtLeftCalorimeter_ARReal_RT_Histo,
488                                k->index, kinEn    488                                k->index, kinEnergy );
489             histoManager->Add( CexmcMissEnFrom    489             histoManager->Add( CexmcMissEnFromLeftCalorimeter_ARReal_RT_Histo,
490                                k->index,          490                                k->index,
491                                kinEnergy - edS    491                                kinEnergy - edStore->calorimeterEDLeft );
492         }                                         492         }
493         if ( tpStore->calorimeterTPRight.IsVal    493         if ( tpStore->calorimeterTPRight.IsValid() )
494         {                                         494         {
495             G4double  kinEnergy( tpStore->calo    495             G4double  kinEnergy( tpStore->calorimeterTPRight.momentumAmp );
496             histoManager->Add( CexmcKinEnAtRig    496             histoManager->Add( CexmcKinEnAtRightCalorimeter_ARReal_RT_Histo,
497                                k->index, kinEn    497                                k->index, kinEnergy );
498             histoManager->Add( CexmcMissEnFrom    498             histoManager->Add( CexmcMissEnFromRightCalorimeter_ARReal_RT_Histo,
499                                k->index,          499                                k->index,
500                                kinEnergy - edS    500                                kinEnergy - edStore->calorimeterEDRight );
501         }                                         501         }
502         if ( tpStore->targetTPOutputParticle.I    502         if ( tpStore->targetTPOutputParticle.IsValid() )
503         {                                         503         {
504             histoManager->Add( CexmcTPInTarget    504             histoManager->Add( CexmcTPInTarget_ARReal_RT_Histo, k->index,
505                            tpStore->targetTPOu    505                            tpStore->targetTPOutputParticle.positionLocal.x(),
506                            tpStore->targetTPOu    506                            tpStore->targetTPOutputParticle.positionLocal.y(),
507                            tpStore->targetTPOu    507                            tpStore->targetTPOutputParticle.positionLocal.z() );
508             histoManager->Add( CexmcKinEnOP_LA    508             histoManager->Add( CexmcKinEnOP_LAB_ARReal_RT_Histo, k->index,
509                                opKinEnergy );     509                                opKinEnergy );
510             histoManager->Add( CexmcAngleOP_SC    510             histoManager->Add( CexmcAngleOP_SCM_ARReal_RT_Histo, k->index,
511                                pmData.outputPa    511                                pmData.outputParticleSCM.cosTheta() );
512             G4double  diffCosTheta( pmData.out    512             G4double  diffCosTheta( pmData.outputParticleSCM.cosTheta() -
513                                     recCosThet    513                                     recCosTheta );
514             histoManager->Add( CexmcDiffAngleO    514             histoManager->Add( CexmcDiffAngleOP_SCM_ARReal_RT_Histo, k->index,
515                                diffCosTheta );    515                                diffCosTheta );
516         }                                         516         }
517         if ( tpStore->targetTPOutputParticleDe    517         if ( tpStore->targetTPOutputParticleDecayProductParticle1.IsValid() &&
518              tpStore->targetTPOutputParticleDe    518              tpStore->targetTPOutputParticleDecayProductParticle2.IsValid() )
519         {                                         519         {
520             G4double  openAngle(                  520             G4double  openAngle(
521                         tpStore->targetTPOutpu    521                         tpStore->targetTPOutputParticleDecayProductParticle1.
522                         directionWorld.angle(     522                         directionWorld.angle( tpStore->
523                                 targetTPOutput    523                                 targetTPOutputParticleDecayProductParticle2.
524                                         direct    524                                         directionWorld ) / deg );
525             histoManager->Add( CexmcOpenAngle_    525             histoManager->Add( CexmcOpenAngle_ARReal_RT_Histo, k->index,
526                                openAngle );       526                                openAngle );
527             G4double  diffOpenAngle( openAngle    527             G4double  diffOpenAngle( openAngle - reconstructor->GetTheAngle() /
528                                      deg );       528                                      deg );
529             histoManager->Add( CexmcDiffOpenAn    529             histoManager->Add( CexmcDiffOpenAngle_ARReal_RT_Histo, k->index,
530                                diffOpenAngle )    530                                diffOpenAngle );
531         }                                         531         }
532         if ( tpStore->calorimeterTPLeft.IsVali    532         if ( tpStore->calorimeterTPLeft.IsValid() )
533         {                                         533         {
534             histoManager->Add( CexmcOPDPAtLeft    534             histoManager->Add( CexmcOPDPAtLeftCalorimeter_ARReal_RT_Histo,
535                                k->index,          535                                k->index,
536                                tpStore->calori    536                                tpStore->calorimeterTPLeft.positionLocal.x(),
537                                tpStore->calori    537                                tpStore->calorimeterTPLeft.positionLocal.y() );
538         }                                         538         }
539         if ( tpStore->calorimeterTPRight.IsVal    539         if ( tpStore->calorimeterTPRight.IsValid() )
540         {                                         540         {
541             histoManager->Add( CexmcOPDPAtRigh    541             histoManager->Add( CexmcOPDPAtRightCalorimeter_ARReal_RT_Histo,
542                                k->index,          542                                k->index,
543                                tpStore->calori    543                                tpStore->calorimeterTPRight.positionLocal.x(),
544                                tpStore->calori    544                                tpStore->calorimeterTPRight.positionLocal.y() );
545         }                                         545         }
546         histoManager->Add( CexmcAbsEnInLeftCal    546         histoManager->Add( CexmcAbsEnInLeftCalorimeter_ARReal_RT_Histo,
547                            k->index, edStore->    547                            k->index, edStore->calorimeterEDLeft );
548         histoManager->Add( CexmcAbsEnInRightCa    548         histoManager->Add( CexmcAbsEnInRightCalorimeter_ARReal_RT_Histo,
549                            k->index, edStore->    549                            k->index, edStore->calorimeterEDRight );
550         histoManager->Add( CexmcRecAngleOP_SCM    550         histoManager->Add( CexmcRecAngleOP_SCM_ARReal_RT_Histo,
551                            k->index, recCosThe    551                            k->index, recCosTheta );
552         histoManager->Add( CexmcRecOpenAngle_A    552         histoManager->Add( CexmcRecOpenAngle_ARReal_RT_Histo,
553                            k->index, reconstru    553                            k->index, reconstructor->GetTheAngle() / deg );
554         histoManager->Add( CexmcRecOPDPAtLeftC    554         histoManager->Add( CexmcRecOPDPAtLeftCalorimeter_ARReal_RT_Histo,
555                            k->index,              555                            k->index,
556                            reconstructor->GetC    556                            reconstructor->GetCalorimeterEPLeftPosition().x(),
557                            reconstructor->GetC    557                            reconstructor->GetCalorimeterEPLeftPosition().y() );
558         histoManager->Add( CexmcRecOPDPAtRight    558         histoManager->Add( CexmcRecOPDPAtRightCalorimeter_ARReal_RT_Histo,
559                            k->index,              559                            k->index,
560                            reconstructor->GetC    560                            reconstructor->GetCalorimeterEPRightPosition().x(),
561                            reconstructor->GetC    561                            reconstructor->GetCalorimeterEPRightPosition().y() );
562     }                                             562     }
563 }                                                 563 }
564                                                   564 
565 #endif                                            565 #endif
566                                                   566 
567                                                   567 
568 void  CexmcEventAction::DrawTrajectories( cons    568 void  CexmcEventAction::DrawTrajectories( const G4Event *  event )
569 {                                                 569 {
570     G4VisManager *  visManager( static_cast< G    570     G4VisManager *  visManager( static_cast< G4VisManager * >(
571                                     G4VVisMana    571                                     G4VVisManager::GetConcreteInstance() ) );
572     if ( ! visManager || ! visManager->GetCurr    572     if ( ! visManager || ! visManager->GetCurrentGraphicsSystem() )
573         return;                                   573         return;
574                                                   574 
575     G4int                    nTraj( 0 );          575     G4int                    nTraj( 0 );
576     G4TrajectoryContainer *  trajContainer( ev    576     G4TrajectoryContainer *  trajContainer( event->GetTrajectoryContainer() );
577                                                   577 
578     if ( ! trajContainer )                        578     if ( ! trajContainer )
579         return;                                   579         return;
580                                                   580 
581     nTraj = trajContainer->entries();             581     nTraj = trajContainer->entries();
582                                                   582 
583     for ( int  i( 0 ); i < nTraj; ++i )           583     for ( int  i( 0 ); i < nTraj; ++i )
584     {                                             584     {
585         G4VTrajectory *  traj( ( *trajContaine    585         G4VTrajectory *  traj( ( *trajContainer )[ i ] );
586         traj->DrawTrajectory();                   586         traj->DrawTrajectory();
587     }                                             587     }
588 }                                                 588 }
589                                                   589 
590                                                   590 
591 void  CexmcEventAction::DrawTrackPoints(          591 void  CexmcEventAction::DrawTrackPoints(
592                                 const CexmcTra    592                                 const CexmcTrackPointsStore *  tpStore ) const
593 {                                                 593 {
594     G4VisManager *  visManager( static_cast< G    594     G4VisManager *  visManager( static_cast< G4VisManager * >(
595                                     G4VVisMana    595                                     G4VVisManager::GetConcreteInstance() ) );
596     if ( ! visManager || ! visManager->GetCurr    596     if ( ! visManager || ! visManager->GetCurrentGraphicsSystem() )
597         return;                                   597         return;
598                                                   598 
599     G4Circle         circle;                      599     G4Circle         circle;
600     G4VisAttributes  visAttributes( CexmcTrack    600     G4VisAttributes  visAttributes( CexmcTrackPointsMarkerColour );
601     circle.SetScreenSize( CexmcSmallCircleScre    601     circle.SetScreenSize( CexmcSmallCircleScreenSize );
602     circle.SetFillStyle( G4Circle::filled );      602     circle.SetFillStyle( G4Circle::filled );
603     circle.SetVisAttributes( visAttributes );     603     circle.SetVisAttributes( visAttributes );
604                                                   604 
605     if ( tpStore->monitorTP.IsValid() )           605     if ( tpStore->monitorTP.IsValid() )
606     {                                             606     {
607         circle.SetPosition( tpStore->monitorTP    607         circle.SetPosition( tpStore->monitorTP.positionWorld );
608         visManager->Draw( circle );               608         visManager->Draw( circle );
609     }                                             609     }
610                                                   610 
611     if ( tpStore->targetTPBeamParticle.IsValid    611     if ( tpStore->targetTPBeamParticle.IsValid() )
612     {                                             612     {
613         circle.SetPosition( tpStore->targetTPB    613         circle.SetPosition( tpStore->targetTPBeamParticle.positionWorld );
614         visManager->Draw( circle );               614         visManager->Draw( circle );
615     }                                             615     }
616                                                   616 
617     if ( tpStore->targetTPOutputParticle.IsVal    617     if ( tpStore->targetTPOutputParticle.IsValid() )
618     {                                             618     {
619         circle.SetPosition( tpStore->targetTPO    619         circle.SetPosition( tpStore->targetTPOutputParticle.positionWorld );
620         visManager->Draw( circle );               620         visManager->Draw( circle );
621     }                                             621     }
622                                                   622 
623     if ( tpStore->vetoCounterTPLeft.IsValid()     623     if ( tpStore->vetoCounterTPLeft.IsValid() )
624     {                                             624     {
625         circle.SetPosition( tpStore->vetoCount    625         circle.SetPosition( tpStore->vetoCounterTPLeft.positionWorld );
626         visManager->Draw( circle );               626         visManager->Draw( circle );
627     }                                             627     }
628                                                   628 
629     if ( tpStore->vetoCounterTPRight.IsValid()    629     if ( tpStore->vetoCounterTPRight.IsValid() )
630     {                                             630     {
631         circle.SetPosition( tpStore->vetoCount    631         circle.SetPosition( tpStore->vetoCounterTPRight.positionWorld );
632         visManager->Draw( circle );               632         visManager->Draw( circle );
633     }                                             633     }
634                                                   634 
635     if ( tpStore->calorimeterTPLeft.IsValid()     635     if ( tpStore->calorimeterTPLeft.IsValid() )
636     {                                             636     {
637         circle.SetPosition( tpStore->calorimet    637         circle.SetPosition( tpStore->calorimeterTPLeft.positionWorld );
638         visManager->Draw( circle );               638         visManager->Draw( circle );
639     }                                             639     }
640                                                   640 
641     if ( tpStore->calorimeterTPRight.IsValid()    641     if ( tpStore->calorimeterTPRight.IsValid() )
642     {                                             642     {
643         circle.SetPosition( tpStore->calorimet    643         circle.SetPosition( tpStore->calorimeterTPRight.positionWorld );
644         visManager->Draw( circle );               644         visManager->Draw( circle );
645     }                                             645     }
646 }                                                 646 }
647                                                   647 
648                                                   648 
649 void  CexmcEventAction::DrawReconstructionData    649 void  CexmcEventAction::DrawReconstructionData( void )
650 {                                                 650 {
651     G4VisManager *  visManager( static_cast< G    651     G4VisManager *  visManager( static_cast< G4VisManager * >(
652                                     G4VVisMana    652                                     G4VVisManager::GetConcreteInstance() ) );
653     if ( ! visManager || ! visManager->GetCurr    653     if ( ! visManager || ! visManager->GetCurrentGraphicsSystem() )
654         return;                                   654         return;
655                                                   655 
656     G4Circle  circle( reconstructor->GetTarget    656     G4Circle  circle( reconstructor->GetTargetEPWorldPosition() );
657     circle.SetScreenSize( CexmcSmallCircleScre    657     circle.SetScreenSize( CexmcSmallCircleScreenSize );
658     circle.SetFillStyle( G4Circle::filled );      658     circle.SetFillStyle( G4Circle::filled );
659     G4VisAttributes  visAttributes( CexmcRecTr    659     G4VisAttributes  visAttributes( CexmcRecTrackPointsMarkerColour );
660     circle.SetVisAttributes( visAttributes );     660     circle.SetVisAttributes( visAttributes );
661     visManager->Draw( circle );                   661     visManager->Draw( circle );
662                                                   662 
663     circle.SetScreenSize( CexmcBigCircleScreen    663     circle.SetScreenSize( CexmcBigCircleScreenSize );
664     circle.SetPosition( reconstructor->GetCalo    664     circle.SetPosition( reconstructor->GetCalorimeterEPLeftWorldPosition() );
665     visManager->Draw( circle );                   665     visManager->Draw( circle );
666                                                   666 
667     circle.SetPosition( reconstructor->GetCalo    667     circle.SetPosition( reconstructor->GetCalorimeterEPRightWorldPosition() );
668     visManager->Draw( circle );                   668     visManager->Draw( circle );
669 }                                                 669 }
670                                                   670 
671                                                   671 
672 void  CexmcEventAction::UpdateRunHits(            672 void  CexmcEventAction::UpdateRunHits(
673                                     const Cexm    673                                     const CexmcAngularRangeList &  aRangesReal,
674                                     const Cexm    674                                     const CexmcAngularRangeList &  aRangesRec,
675                                     G4bool  tp    675                                     G4bool  tpDigitizerHasTriggered,
676                                     G4bool  ed    676                                     G4bool  edDigitizerHasTriggered,
677                                     G4bool  ed    677                                     G4bool  edDigitizerMonitorHasTriggered,
678                                     G4bool  re    678                                     G4bool  reconstructorHasFullTrigger,
679                                     const Cexm    679                                     const CexmcAngularRange &  aGap )
680 {                                                 680 {
681     G4RunManager *    runManager( G4RunManager    681     G4RunManager *    runManager( G4RunManager::GetRunManager() );
682     const CexmcRun *  run( static_cast< const     682     const CexmcRun *  run( static_cast< const CexmcRun * >(
683                                                   683                                                 runManager->GetCurrentRun() ) );
684     CexmcRun *        theRun( const_cast< Cexm    684     CexmcRun *        theRun( const_cast< CexmcRun * >( run ) );
685                                                   685 
686     if ( tpDigitizerHasTriggered )                686     if ( tpDigitizerHasTriggered )
687     {                                             687     {
688         for ( CexmcAngularRangeList::const_ite    688         for ( CexmcAngularRangeList::const_iterator  k( aRangesReal.begin() );
689                                                   689                                                 k != aRangesReal.end(); ++k )
690         {                                         690         {
691             theRun->IncrementNmbOfHitsSampledF    691             theRun->IncrementNmbOfHitsSampledFull( k->index );
692             if ( edDigitizerMonitorHasTriggere    692             if ( edDigitizerMonitorHasTriggered )
693                 theRun->IncrementNmbOfHitsSamp    693                 theRun->IncrementNmbOfHitsSampled( k->index );
694             if ( reconstructorHasFullTrigger )    694             if ( reconstructorHasFullTrigger )
695                 theRun->IncrementNmbOfHitsTrig    695                 theRun->IncrementNmbOfHitsTriggeredRealRange( k->index );
696         }                                         696         }
697         if ( reconstructorHasFullTrigger )        697         if ( reconstructorHasFullTrigger )
698         {                                         698         {
699             if ( aRangesRec.empty() )             699             if ( aRangesRec.empty() )
700             {                                     700             {
701                 theRun->IncrementNmbOfOrphanHi    701                 theRun->IncrementNmbOfOrphanHits( aGap.index );
702             }                                     702             }
703             else                                  703             else
704             {                                     704             {
705                 for ( CexmcAngularRangeList::c    705                 for ( CexmcAngularRangeList::const_iterator
706                         k( aRangesRec.begin()     706                         k( aRangesRec.begin() ); k != aRangesRec.end(); ++k )
707                 {                                 707                 {
708                     theRun->IncrementNmbOfHits    708                     theRun->IncrementNmbOfHitsTriggeredRecRange( k->index );
709                 }                                 709                 }
710             }                                     710             }
711         }                                         711         }
712     }                                             712     }
713     else                                          713     else
714     {                                             714     {
715         if ( edDigitizerHasTriggered )            715         if ( edDigitizerHasTriggered )
716             theRun->IncrementNmbOfFalseHitsTri    716             theRun->IncrementNmbOfFalseHitsTriggeredEDT();
717         if ( reconstructorHasFullTrigger )        717         if ( reconstructorHasFullTrigger )
718             theRun->IncrementNmbOfFalseHitsTri    718             theRun->IncrementNmbOfFalseHitsTriggeredRec();
719     }                                             719     }
720 }                                                 720 }
721                                                   721 
722                                                   722 
723 #ifdef CEXMC_USE_PERSISTENCY                      723 #ifdef CEXMC_USE_PERSISTENCY
724                                                   724 
725 void  CexmcEventAction::SaveEvent( const G4Eve    725 void  CexmcEventAction::SaveEvent( const G4Event *  event,
726                                    G4bool  edD    726                                    G4bool  edDigitizerHasTriggered,
727                                    const Cexmc    727                                    const CexmcEnergyDepositStore *  edStore,
728                                    const Cexmc    728                                    const CexmcTrackPointsStore *  tpStore,
729                                    const Cexmc    729                                    const CexmcProductionModelData &  pmData )
730 {                                                 730 {
731     CexmcRunManager *  runManager( static_cast    731     CexmcRunManager *  runManager( static_cast< CexmcRunManager * >(
732                                             G4    732                                             G4RunManager::GetRunManager() ) );
733     if ( ! runManager->ProjectIsSaved() )         733     if ( ! runManager->ProjectIsSaved() )
734         return;                                   734         return;
735                                                   735 
736     if ( runManager->GetEventDataVerboseLevel(    736     if ( runManager->GetEventDataVerboseLevel() == CexmcWriteNoEventData )
737         return;                                   737         return;
738                                                   738 
739     if ( ! edDigitizerHasTriggered && runManag    739     if ( ! edDigitizerHasTriggered && runManager->GetEventDataVerboseLevel() !=
740                                                   740                                                 CexmcWriteEventDataOnEveryTPT )
741         return;                                   741         return;
742                                                   742 
743     boost::archive::binary_oarchive *  archive    743     boost::archive::binary_oarchive *  archive(
744                                             ru    744                                             runManager->GetEventsArchive() );
745     if ( archive )                                745     if ( archive )
746     {                                             746     {
747         CexmcEventSObject  sObject = { event->    747         CexmcEventSObject  sObject = { event->GetEventID(),
748             edDigitizerHasTriggered, edStore->    748             edDigitizerHasTriggered, edStore->monitorED,
749             edStore->vetoCounterEDLeft, edStor    749             edStore->vetoCounterEDLeft, edStore->vetoCounterEDRight,
750             edStore->calorimeterEDLeft, edStor    750             edStore->calorimeterEDLeft, edStore->calorimeterEDRight,
751             edStore->calorimeterEDLeftCollecti    751             edStore->calorimeterEDLeftCollection,
752             edStore->calorimeterEDRightCollect    752             edStore->calorimeterEDRightCollection,
753             tpStore->monitorTP, tpStore->targe    753             tpStore->monitorTP, tpStore->targetTPBeamParticle,
754             tpStore->targetTPOutputParticle, t    754             tpStore->targetTPOutputParticle, tpStore->targetTPNucleusParticle,
755             tpStore->targetTPOutputParticleDec    755             tpStore->targetTPOutputParticleDecayProductParticle1,
756             tpStore->targetTPOutputParticleDec    756             tpStore->targetTPOutputParticleDecayProductParticle2,
757             tpStore->vetoCounterTPLeft, tpStor    757             tpStore->vetoCounterTPLeft, tpStore->vetoCounterTPRight,
758             tpStore->calorimeterTPLeft, tpStor    758             tpStore->calorimeterTPLeft, tpStore->calorimeterTPRight, pmData };
759         archive->operator<<( sObject );           759         archive->operator<<( sObject );
760         const CexmcRun *  run( static_cast< co    760         const CexmcRun *  run( static_cast< const CexmcRun * >(
761                                                   761                                                 runManager->GetCurrentRun() ) );
762         CexmcRun *        theRun( const_cast<     762         CexmcRun *        theRun( const_cast< CexmcRun * >( run ) );
763         theRun->IncrementNmbOfSavedEvents();      763         theRun->IncrementNmbOfSavedEvents();
764     }                                             764     }
765 }                                                 765 }
766                                                   766 
767                                                   767 
768 void  CexmcEventAction::SaveEventFast( const G    768 void  CexmcEventAction::SaveEventFast( const G4Event *  event,
769                                        G4bool     769                                        G4bool  tpDigitizerHasTriggered,
770                                        G4bool     770                                        G4bool  edDigitizerHasTriggered,
771                                        G4bool     771                                        G4bool  edDigitizerMonitorHasTriggered,
772                                        G4doubl    772                                        G4double  opCosThetaSCM )
773 {                                                 773 {
774     CexmcRunManager *  runManager( static_cast    774     CexmcRunManager *  runManager( static_cast< CexmcRunManager * >(
775                                             G4    775                                             G4RunManager::GetRunManager() ) );
776     if ( ! runManager->ProjectIsSaved() )         776     if ( ! runManager->ProjectIsSaved() )
777         return;                                   777         return;
778                                                   778 
779     if ( runManager->GetEventDataVerboseLevel(    779     if ( runManager->GetEventDataVerboseLevel() == CexmcWriteNoEventData )
780         return;                                   780         return;
781                                                   781 
782     boost::archive::binary_oarchive *  archive    782     boost::archive::binary_oarchive *  archive(
783                                         runMan    783                                         runManager->GetFastEventsArchive() );
784     if ( archive )                                784     if ( archive )
785     {                                             785     {
786         if ( ! tpDigitizerHasTriggered )          786         if ( ! tpDigitizerHasTriggered )
787             opCosThetaSCM = CexmcInvalidCosThe    787             opCosThetaSCM = CexmcInvalidCosTheta;
788                                                   788 
789         CexmcEventFastSObject  sObject = { eve    789         CexmcEventFastSObject  sObject = { event->GetEventID(), opCosThetaSCM,
790                                            edD    790                                            edDigitizerHasTriggered,
791                                            edD    791                                            edDigitizerMonitorHasTriggered };
792         archive->operator<<( sObject );           792         archive->operator<<( sObject );
793         const CexmcRun *  run( static_cast< co    793         const CexmcRun *  run( static_cast< const CexmcRun * >(
794                                                   794                                                 runManager->GetCurrentRun() ) );
795         CexmcRun *        theRun( const_cast<     795         CexmcRun *        theRun( const_cast< CexmcRun * >( run ) );
796         theRun->IncrementNmbOfSavedFastEvents(    796         theRun->IncrementNmbOfSavedFastEvents();
797     }                                             797     }
798 }                                                 798 }
799                                                   799 
800 #endif                                            800 #endif
801                                                   801 
802                                                   802 
803 void  CexmcEventAction::EndOfEventAction( cons    803 void  CexmcEventAction::EndOfEventAction( const G4Event *  event )
804 {                                                 804 {
805     G4DigiManager *                digiManager    805     G4DigiManager *                digiManager( G4DigiManager::GetDMpointer() );
806     CexmcEnergyDepositDigitizer *  energyDepos    806     CexmcEnergyDepositDigitizer *  energyDepositDigitizer(
807             static_cast< CexmcEnergyDepositDig    807             static_cast< CexmcEnergyDepositDigitizer * >( digiManager->
808                                 FindDigitizerM    808                                 FindDigitizerModule( CexmcEDDigitizerName ) ) );
809     CexmcTrackPointsDigitizer *    trackPoints    809     CexmcTrackPointsDigitizer *    trackPointsDigitizer(
810             static_cast< CexmcTrackPointsDigit    810             static_cast< CexmcTrackPointsDigitizer * >( digiManager->
811                                 FindDigitizerM    811                                 FindDigitizerModule( CexmcTPDigitizerName ) ) );
812                                                   812 
813     energyDepositDigitizer->Digitize();           813     energyDepositDigitizer->Digitize();
814     trackPointsDigitizer->Digitize();             814     trackPointsDigitizer->Digitize();
815                                                   815 
816     G4bool  edDigitizerMonitorHasTriggered(       816     G4bool  edDigitizerMonitorHasTriggered(
817                                 energyDepositD    817                                 energyDepositDigitizer->MonitorHasTriggered() );
818     G4bool  edDigitizerHasTriggered( false );     818     G4bool  edDigitizerHasTriggered( false );
819                                                   819 
820     CexmcEventInfo *  eventInfo( static_cast<     820     CexmcEventInfo *  eventInfo( static_cast< CexmcEventInfo * >(
821                                                   821                                                 event->GetUserInformation() ) );
822     if ( ! eventInfo || eventInfo->EdTriggerIs    822     if ( ! eventInfo || eventInfo->EdTriggerIsOk() )
823         edDigitizerHasTriggered = energyDeposi    823         edDigitizerHasTriggered = energyDepositDigitizer->HasTriggered();
824                                                   824 
825     G4bool  tpDigitizerHasTriggered( trackPoin    825     G4bool  tpDigitizerHasTriggered( trackPointsDigitizer->HasTriggered() );
826     G4bool  reconstructorHasBasicTrigger( fals    826     G4bool  reconstructorHasBasicTrigger( false );
827     G4bool  reconstructorHasFullTrigger( false    827     G4bool  reconstructorHasFullTrigger( false );
828                                                   828 
829     CexmcEnergyDepositStore *  edStore( MakeEn    829     CexmcEnergyDepositStore *  edStore( MakeEnergyDepositStore(
830                                                   830                                                     energyDepositDigitizer ) );
831     CexmcTrackPointsStore *    tpStore( MakeTr    831     CexmcTrackPointsStore *    tpStore( MakeTrackPointsStore(
832                                                   832                                                     trackPointsDigitizer ) );
833                                                   833 
834     try                                           834     try
835     {                                             835     {
836         CexmcProductionModel *  productionMode    836         CexmcProductionModel *  productionModel(
837                                         physic    837                                         physicsManager->GetProductionModel() );
838                                                   838 
839         if ( ! productionModel )                  839         if ( ! productionModel )
840             throw CexmcException( CexmcWeirdEx    840             throw CexmcException( CexmcWeirdException );
841                                                   841 
842         const CexmcAngularRangeList &     angu    842         const CexmcAngularRangeList &     angularRanges(
843                                 productionMode    843                                 productionModel->GetAngularRanges() );
844         const CexmcAngularRangeList &     trig    844         const CexmcAngularRangeList &     triggeredAngularRanges(
845                                 productionMode    845                                 productionModel->GetTriggeredAngularRanges() );
846         const CexmcProductionModelData &  pmDa    846         const CexmcProductionModelData &  pmData(
847                                 productionMode    847                                 productionModel->GetProductionModelData() );
848                                                   848 
849         if ( edDigitizerHasTriggered )            849         if ( edDigitizerHasTriggered )
850         {                                         850         {
851             reconstructor->Reconstruct( edStor    851             reconstructor->Reconstruct( edStore );
852             reconstructorHasBasicTrigger = rec    852             reconstructorHasBasicTrigger = reconstructor->HasBasicTrigger();
853             reconstructorHasFullTrigger = reco    853             reconstructorHasFullTrigger = reconstructor->HasFullTrigger();
854         }                                         854         }
855                                                   855 
856         CexmcAngularRangeList  triggeredRecAng    856         CexmcAngularRangeList  triggeredRecAngularRanges;
857                                                   857 
858         if ( reconstructorHasBasicTrigger )       858         if ( reconstructorHasBasicTrigger )
859         {                                         859         {
860             for ( CexmcAngularRangeList::const    860             for ( CexmcAngularRangeList::const_iterator
861                   k( angularRanges.begin() );     861                   k( angularRanges.begin() ); k != angularRanges.end(); ++k )
862             {                                     862             {
863                 G4double  cosTheta( reconstruc    863                 G4double  cosTheta( reconstructor->GetProductionModelData().
864                                     outputPart    864                                     outputParticleSCM.cosTheta() );
865                 if ( cosTheta <= k->top && cos    865                 if ( cosTheta <= k->top && cosTheta > k->bottom )
866                     triggeredRecAngularRanges.    866                     triggeredRecAngularRanges.push_back( CexmcAngularRange(
867                                                   867                                                 k->top, k->bottom, k->index ) );
868             }                                     868             }
869         }                                         869         }
870                                                   870 
871         CexmcAngularRange  angularGap( 0.0, 0.    871         CexmcAngularRange  angularGap( 0.0, 0.0, 0 );
872         if ( triggeredRecAngularRanges.empty()    872         if ( triggeredRecAngularRanges.empty() )
873         {                                         873         {
874             CexmcAngularRangeList  angularGaps    874             CexmcAngularRangeList  angularGaps;
875             GetAngularGaps( angularRanges, ang    875             GetAngularGaps( angularRanges, angularGaps );
876             for ( CexmcAngularRangeList::const    876             for ( CexmcAngularRangeList::const_iterator
877                     k( angularGaps.begin() );     877                     k( angularGaps.begin() ); k != angularGaps.end(); ++k )
878             {                                     878             {
879                 G4double  cosTheta( reconstruc    879                 G4double  cosTheta( reconstructor->GetProductionModelData().
880                                     outputPart    880                                     outputParticleSCM.cosTheta() );
881                 if ( cosTheta <= k->top && cos    881                 if ( cosTheta <= k->top && cosTheta > k->bottom )
882                 {                                 882                 {
883                     angularGap = *k;              883                     angularGap = *k;
884                     break;                        884                     break;
885                 }                                 885                 }
886             }                                     886             }
887         }                                         887         }
888                                                   888 
889         UpdateRunHits( triggeredAngularRanges,    889         UpdateRunHits( triggeredAngularRanges, triggeredRecAngularRanges,
890                        tpDigitizerHasTriggered    890                        tpDigitizerHasTriggered, edDigitizerHasTriggered,
891                        edDigitizerMonitorHasTr    891                        edDigitizerMonitorHasTriggered,
892                        reconstructorHasFullTri    892                        reconstructorHasFullTrigger, angularGap );
893                                                   893 
894         if ( verbose > 0 )                        894         if ( verbose > 0 )
895         {                                         895         {
896             G4bool  printMessages( verbose > 3    896             G4bool  printMessages( verbose > 3 ||
897                         ( ( verbose == 1 ) &&     897                         ( ( verbose == 1 ) && tpDigitizerHasTriggered ) ||
898                         ( ( verbose == 2 ) &&     898                         ( ( verbose == 2 ) && edDigitizerHasTriggered ) ||
899                         ( ( verbose == 3 ) &&     899                         ( ( verbose == 3 ) && ( tpDigitizerHasTriggered ||
900                                                   900                                                 edDigitizerHasTriggered ) ) );
901             if ( printMessages )                  901             if ( printMessages )
902             {                                     902             {
903                 G4cout << "Event " << event->G    903                 G4cout << "Event " << event->GetEventID() << G4endl;
904                 if ( tpDigitizerHasTriggered )    904                 if ( tpDigitizerHasTriggered )
905                 {                                 905                 {
906                     PrintTrackPoints( tpStore     906                     PrintTrackPoints( tpStore );
907                     PrintProductionModelData(     907                     PrintProductionModelData( triggeredAngularRanges, pmData );
908                 }                                 908                 }
909                 if ( reconstructorHasBasicTrig    909                 if ( reconstructorHasBasicTrigger )
910                     PrintReconstructedData( tr    910                     PrintReconstructedData( triggeredRecAngularRanges,
911                                             an    911                                             angularGap );
912                 if ( edDigitizerHasTriggered )    912                 if ( edDigitizerHasTriggered )
913                     PrintEnergyDeposit( edStor    913                     PrintEnergyDeposit( edStore );
914             }                                     914             }
915         }                                         915         }
916                                                   916 
917         if ( verboseDraw > 0 )                    917         if ( verboseDraw > 0 )
918         {                                         918         {
919             G4bool  drawTrajectories( verboseD    919             G4bool  drawTrajectories( verboseDraw > 3 ||
920                         ( ( verboseDraw == 1 )    920                         ( ( verboseDraw == 1 ) && tpDigitizerHasTriggered ) ||
921                         ( ( verboseDraw == 2 )    921                         ( ( verboseDraw == 2 ) && edDigitizerHasTriggered ) ||
922                         ( ( verboseDraw == 3 )    922                         ( ( verboseDraw == 3 ) && ( tpDigitizerHasTriggered ||
923                                                   923                                                 edDigitizerHasTriggered ) ) );
924             if ( drawTrajectories )               924             if ( drawTrajectories )
925             {                                     925             {
926                 DrawTrajectories( event );        926                 DrawTrajectories( event );
927                 if ( tpDigitizerHasTriggered )    927                 if ( tpDigitizerHasTriggered )
928                     DrawTrackPoints( tpStore )    928                     DrawTrackPoints( tpStore );
929                 if ( reconstructorHasBasicTrig    929                 if ( reconstructorHasBasicTrigger )
930                     DrawReconstructionData();     930                     DrawReconstructionData();
931             }                                     931             }
932         }                                         932         }
933                                                   933 
934 #ifdef CEXMC_USE_PERSISTENCY                      934 #ifdef CEXMC_USE_PERSISTENCY
935         if ( edDigitizerHasTriggered || tpDigi    935         if ( edDigitizerHasTriggered || tpDigitizerHasTriggered )
936         {                                         936         {
937             SaveEventFast( event, tpDigitizerH    937             SaveEventFast( event, tpDigitizerHasTriggered,
938                            edDigitizerHasTrigg    938                            edDigitizerHasTriggered,
939                            edDigitizerMonitorH    939                            edDigitizerMonitorHasTriggered,
940                            pmData.outputPartic    940                            pmData.outputParticleSCM.cosTheta() );
941             SaveEvent( event, edDigitizerHasTr    941             SaveEvent( event, edDigitizerHasTriggered, edStore, tpStore,
942                        pmData );                  942                        pmData );
943         }                                         943         }
944 #endif                                            944 #endif
945                                                   945 
946 #ifdef CEXMC_USE_ROOT                             946 #ifdef CEXMC_USE_ROOT
947         /* opKinEnergy will be used in several    947         /* opKinEnergy will be used in several histos */
948         if ( tpStore->targetTPOutputParticle.I    948         if ( tpStore->targetTPOutputParticle.IsValid() )
949         {                                         949         {
950             opKinEnergy = CexmcGetKinEnergy(      950             opKinEnergy = CexmcGetKinEnergy(
951                     tpStore->targetTPOutputPar    951                     tpStore->targetTPOutputParticle.momentumAmp,
952                     tpStore->targetTPOutputPar    952                     tpStore->targetTPOutputParticle.particle->GetPDGMass() );
953         }                                         953         }
954                                                   954 
955         if ( edDigitizerHasTriggered )            955         if ( edDigitizerHasTriggered )
956             FillEDTHistos( edStore, triggeredA    956             FillEDTHistos( edStore, triggeredAngularRanges );
957                                                   957 
958         /* fill TPT histos only when the monit    958         /* fill TPT histos only when the monitor has triggered because events
959          * when it was missed have less value     959          * when it was missed have less value for us */
960         if ( tpDigitizerHasTriggered && edDigi    960         if ( tpDigitizerHasTriggered && edDigitizerMonitorHasTriggered )
961             FillTPTHistos( tpStore, pmData, tr    961             FillTPTHistos( tpStore, pmData, triggeredAngularRanges );
962                                                   962 
963         if ( reconstructorHasBasicTrigger )       963         if ( reconstructorHasBasicTrigger )
964             FillRTHistos( reconstructorHasFull    964             FillRTHistos( reconstructorHasFullTrigger, edStore, tpStore,
965                           pmData, triggeredAng    965                           pmData, triggeredAngularRanges );
966 #endif                                            966 #endif
967                                                   967 
968         G4Event *  theEvent( const_cast< G4Eve    968         G4Event *  theEvent( const_cast< G4Event * >( event ) );
969         if ( eventInfo )                          969         if ( eventInfo )
970         {                                         970         {
971             delete eventInfo;                     971             delete eventInfo;
972             theEvent->SetUserInformation( NULL    972             theEvent->SetUserInformation( NULL );
973         }                                         973         }
974         theEvent->SetUserInformation( new Cexm    974         theEvent->SetUserInformation( new CexmcEventInfo(
975                                                   975                                                 edDigitizerHasTriggered,
976                                                   976                                                 tpDigitizerHasTriggered,
977                                                   977                                                 reconstructorHasFullTrigger ) );
978     }                                             978     }
979     catch ( CexmcException &  e )                 979     catch ( CexmcException &  e )
980     {                                             980     {
981         G4cout << e.what() << G4endl;             981         G4cout << e.what() << G4endl;
982     }                                             982     }
983     catch ( ... )                                 983     catch ( ... )
984     {                                             984     {
985         G4cout << "Unknown exception caught" <    985         G4cout << "Unknown exception caught" << G4endl;
986     }                                             986     }
987                                                   987 
988     delete edStore;                               988     delete edStore;
989     delete tpStore;                               989     delete tpStore;
990 }                                                 990 }
991                                                   991 
992                                                   992