| Geant4 Cross Reference |
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: CexmcHistoManager.cc 29 * Filename: CexmcHistoManager.cc
30 * 30 *
31 * Description: histograming manager (sing 31 * Description: histograming manager (singleton)
32 * 32 *
33 * Version: 1.0 33 * Version: 1.0
34 * Created: 26.11.2009 21:00:03 34 * Created: 26.11.2009 21:00:03
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 #ifdef CEXMC_USE_ROOT 44 #ifdef CEXMC_USE_ROOT
45 45
46 #include <iostream> 46 #include <iostream>
47 #include <iomanip> 47 #include <iomanip>
48 #include <TH1.h> 48 #include <TH1.h>
49 #include <TH1F.h> 49 #include <TH1F.h>
50 #include <TH2F.h> 50 #include <TH2F.h>
51 #include <TH3F.h> 51 #include <TH3F.h>
52 #include <TFile.h> 52 #include <TFile.h>
53 #include <TObject.h> 53 #include <TObject.h>
54 #include <TCollection.h> 54 #include <TCollection.h>
55 #include <TDirectory.h> 55 #include <TDirectory.h>
56 #include <TString.h> 56 #include <TString.h>
57 #ifdef CEXMC_USE_ROOTQT 57 #ifdef CEXMC_USE_ROOTQT
58 #include <TCanvas.h> 58 #include <TCanvas.h>
59 #include <TList.h> 59 #include <TList.h>
60 #include <QApplication> 60 #include <QApplication>
61 #include <QFont> 61 #include <QFont>
62 #include <QMenu> 62 #include <QMenu>
63 #include <G4UIQt.hh> 63 #include <G4UIQt.hh>
64 #include "CexmcMessenger.hh" 64 #include "CexmcMessenger.hh"
65 #endif 65 #endif
66 #include <G4LogicalVolume.hh> 66 #include <G4LogicalVolume.hh>
67 #include <G4Box.hh> 67 #include <G4Box.hh>
68 #include <G4Tubs.hh> 68 #include <G4Tubs.hh>
69 #include <G4SystemOfUnits.hh> 69 #include <G4SystemOfUnits.hh>
70 #include "CexmcHistoManager.hh" 70 #include "CexmcHistoManager.hh"
71 #include "CexmcHistoManagerMessenger.hh" 71 #include "CexmcHistoManagerMessenger.hh"
72 #include "CexmcProductionModel.hh" 72 #include "CexmcProductionModel.hh"
73 #include "CexmcPhysicsManager.hh" 73 #include "CexmcPhysicsManager.hh"
74 #include "CexmcRunManager.hh" 74 #include "CexmcRunManager.hh"
75 #include "CexmcSetup.hh" 75 #include "CexmcSetup.hh"
76 #include "CexmcException.hh" 76 #include "CexmcException.hh"
77 #include "CexmcHistoWidget.hh" 77 #include "CexmcHistoWidget.hh"
78 78
79 79
80 namespace 80 namespace
81 { 81 {
82 const G4double CexmcHistoBeamMomentumMin( 82 const G4double CexmcHistoBeamMomentumMin( 0.0 * GeV );
83 const G4double CexmcHistoBeamMomentumMax( 83 const G4double CexmcHistoBeamMomentumMax( 1.0 * GeV );
84 const G4double CexmcHistoBeamMomentumReso 84 const G4double CexmcHistoBeamMomentumResolution( 0.5 * MeV );
85 const G4double CexmcHistoTPResolution( 0. 85 const G4double CexmcHistoTPResolution( 0.1 * cm );
86 const G4double CexmcHistoTPSafetyArea( 1. 86 const G4double CexmcHistoTPSafetyArea( 1.0 * cm );
87 const G4double CexmcHistoMassResolution( 87 const G4double CexmcHistoMassResolution( 1.0 * MeV );
88 const G4double CexmcHistoEnergyMax( 1.0 * 88 const G4double CexmcHistoEnergyMax( 1.0 * GeV );
89 const G4double CexmcHistoEnergyResolution 89 const G4double CexmcHistoEnergyResolution( 1.0 * MeV );
90 const G4double CexmcHistoMissEnergyMin( - 90 const G4double CexmcHistoMissEnergyMin( -0.1 * GeV );
91 const G4double CexmcHistoMissEnergyMax( 0 91 const G4double CexmcHistoMissEnergyMax( 0.2 * GeV );
92 const G4double CexmcHistoMissEnergyResolu 92 const G4double CexmcHistoMissEnergyResolution( 0.2 * MeV );
93 const G4double CexmcHistoAngularResolutio 93 const G4double CexmcHistoAngularResolution( 0.5 );
94 const G4double CexmcHistoAngularCResoluti 94 const G4double CexmcHistoAngularCResolution( 0.001 );
95 #ifdef CEXMC_USE_ROOTQT 95 #ifdef CEXMC_USE_ROOTQT
96 const G4int CexmcHistoCanvasWidth( 800 96 const G4int CexmcHistoCanvasWidth( 800 );
97 const G4int CexmcHistoCanvasHeight( 60 97 const G4int CexmcHistoCanvasHeight( 600 );
98 #endif 98 #endif
99 const G4String CexmcHistoDirectoryHandle( 99 const G4String CexmcHistoDirectoryHandle( "histograms" );
100 const G4String CexmcHistoDirectoryTitle( 100 const G4String CexmcHistoDirectoryTitle( "Histograms" );
101 } 101 }
102 102
103 103
104 CexmcHistoManager * CexmcHistoManager::instan 104 CexmcHistoManager * CexmcHistoManager::instance( NULL );
105 105
106 106
107 CexmcHistoManager * CexmcHistoManager::Instan 107 CexmcHistoManager * CexmcHistoManager::Instance( void )
108 { 108 {
109 if ( instance == NULL ) 109 if ( instance == NULL )
110 instance = new CexmcHistoManager; 110 instance = new CexmcHistoManager;
111 111
112 return instance; 112 return instance;
113 } 113 }
114 114
115 115
116 void CexmcHistoManager::Destroy( void ) 116 void CexmcHistoManager::Destroy( void )
117 { 117 {
118 delete instance; 118 delete instance;
119 instance = NULL; 119 instance = NULL;
120 } 120 }
121 121
122 122
123 CexmcHistoManager::CexmcHistoManager() : outFi 123 CexmcHistoManager::CexmcHistoManager() : outFile( NULL ),
124 isInitialized( false ), opName( "" ), nopN 124 isInitialized( false ), opName( "" ), nopName( "" ), opMass( 0. ),
125 nopMass( 0. ), verboseLevel( 0 ), 125 nopMass( 0. ), verboseLevel( 0 ),
126 #ifdef CEXMC_USE_ROOTQT 126 #ifdef CEXMC_USE_ROOTQT
127 rootCanvas( NULL ), areLiveHistogramsEnabl 127 rootCanvas( NULL ), areLiveHistogramsEnabled( false ),
128 isHistoMenuInitialized( false ), drawOptio 128 isHistoMenuInitialized( false ), drawOptions1D( "" ), drawOptions2D( "" ),
129 drawOptions3D( "" ), histoMenuHandle( "" ) 129 drawOptions3D( "" ), histoMenuHandle( "" ), histoMenuLabel( "" ),
130 #endif 130 #endif
131 messenger( NULL ) 131 messenger( NULL )
132 { 132 {
133 for ( int i( 0 ); i < CexmcHistoType_SIZE 133 for ( int i( 0 ); i < CexmcHistoType_SIZE; ++i )
134 { 134 {
135 histos.insert( CexmcHistoPair( CexmcHi 135 histos.insert( CexmcHistoPair( CexmcHistoType( i ),
136 CexmcHi 136 CexmcHistoVector() ) );
137 } 137 }
138 138
139 messenger = new CexmcHistoManagerMessenger 139 messenger = new CexmcHistoManagerMessenger( this );
140 } 140 }
141 141
142 142
143 CexmcHistoManager::~CexmcHistoManager() 143 CexmcHistoManager::~CexmcHistoManager()
144 { 144 {
145 if ( outFile ) 145 if ( outFile )
146 { 146 {
147 outFile->Write(); 147 outFile->Write();
148 outFile->Close(); 148 outFile->Close();
149 } 149 }
150 150
151 /* all histograms will be deleted by outFi 151 /* all histograms will be deleted by outFile destructor! */
152 delete outFile; 152 delete outFile;
153 #ifdef CEXMC_USE_ROOTQT 153 #ifdef CEXMC_USE_ROOTQT
154 delete rootCanvas; 154 delete rootCanvas;
155 #endif 155 #endif
156 delete messenger; 156 delete messenger;
157 } 157 }
158 158
159 159
160 void CexmcHistoManager::AddHisto( const Cexmc 160 void CexmcHistoManager::AddHisto( const CexmcHistoData & data,
161 const Cexmc 161 const CexmcAngularRange & aRange )
162 { 162 {
163 G4String fullName( data.name ); 163 G4String fullName( data.name );
164 G4String fullTitle( data.title ); 164 G4String fullTitle( data.title );
165 G4String rangeTypeLabel; 165 G4String rangeTypeLabel;
166 G4String triggerTypeLabel; 166 G4String triggerTypeLabel;
167 G4String decorTriggerTypeLabel; 167 G4String decorTriggerTypeLabel;
168 168
169 if ( data.isARHisto ) 169 if ( data.isARHisto )
170 { 170 {
171 if ( data.isARRec ) 171 if ( data.isARRec )
172 { 172 {
173 fullName += "_arrec"; 173 fullName += "_arrec";
174 rangeTypeLabel = "rec"; 174 rangeTypeLabel = "rec";
175 } 175 }
176 else 176 else
177 { 177 {
178 fullName += "_arreal"; 178 fullName += "_arreal";
179 rangeTypeLabel = "real"; 179 rangeTypeLabel = "real";
180 } 180 }
181 } 181 }
182 182
183 switch ( data.triggerType ) 183 switch ( data.triggerType )
184 { 184 {
185 case CexmcTPT : 185 case CexmcTPT :
186 fullName += "_tpt"; 186 fullName += "_tpt";
187 decorTriggerTypeLabel = " --tpt--"; 187 decorTriggerTypeLabel = " --tpt--";
188 fullTitle += decorTriggerTypeLabel; 188 fullTitle += decorTriggerTypeLabel;
189 triggerTypeLabel = "tpt"; 189 triggerTypeLabel = "tpt";
190 break; 190 break;
191 case CexmcEDT : 191 case CexmcEDT :
192 fullName += "_edt"; 192 fullName += "_edt";
193 decorTriggerTypeLabel = " --edt--"; 193 decorTriggerTypeLabel = " --edt--";
194 fullTitle += decorTriggerTypeLabel; 194 fullTitle += decorTriggerTypeLabel;
195 triggerTypeLabel = "edt"; 195 triggerTypeLabel = "edt";
196 break; 196 break;
197 case CexmcRT : 197 case CexmcRT :
198 fullName += "_rt"; 198 fullName += "_rt";
199 decorTriggerTypeLabel = " --rt--"; 199 decorTriggerTypeLabel = " --rt--";
200 fullTitle += decorTriggerTypeLabel; 200 fullTitle += decorTriggerTypeLabel;
201 triggerTypeLabel = "rt"; 201 triggerTypeLabel = "rt";
202 break; 202 break;
203 default : 203 default :
204 break; 204 break;
205 } 205 }
206 206
207 CexmcHistosMap::iterator found( histos.fi 207 CexmcHistosMap::iterator found( histos.find( data.type ) );
208 208
209 if ( found == histos.end() ) 209 if ( found == histos.end() )
210 throw CexmcException( CexmcWeirdExcept 210 throw CexmcException( CexmcWeirdException );
211 211
212 CexmcHistoVector & histoVector( found->se 212 CexmcHistoVector & histoVector( found->second );
213 213
214 if ( data.isARHisto ) 214 if ( data.isARHisto )
215 { 215 {
216 G4bool dirOk( false ); 216 G4bool dirOk( false );
217 217
218 dirOk = gDirectory->Get( fullName ) != 218 dirOk = gDirectory->Get( fullName ) != NULL;
219 219
220 if ( ! dirOk ) 220 if ( ! dirOk )
221 dirOk = ( gDirectory->mkdir( fullN 221 dirOk = ( gDirectory->mkdir( fullName, fullTitle ) != NULL );
222 222
223 if ( dirOk ) 223 if ( dirOk )
224 gDirectory->cd( fullName ); 224 gDirectory->cd( fullName );
225 225
226 std::ostringstream histoName; 226 std::ostringstream histoName;
227 std::ostringstream histoTitle; 227 std::ostringstream histoTitle;
228 histoName << data.name << "_r" << aRan 228 histoName << data.name << "_r" << aRange.index + 1 << rangeTypeLabel <<
229 "_" << triggerTypeLabel; 229 "_" << triggerTypeLabel;
230 histoTitle << data.title << " {range " 230 histoTitle << data.title << " {range " << aRange.index + 1 <<
231 rangeTypeLabel << " [" << std: 231 rangeTypeLabel << " [" << std::fixed <<
232 std::setprecision( 4 ) << aRan 232 std::setprecision( 4 ) << aRange.top << ", " <<
233 aRange.bottom << ")}" << decor 233 aRange.bottom << ")}" << decorTriggerTypeLabel;
234 CreateHisto( histoVector, data.impl, h 234 CreateHisto( histoVector, data.impl, histoName.str(), histoTitle.str(),
235 data.axes ); 235 data.axes );
236 236
237 if ( dirOk ) 237 if ( dirOk )
238 gDirectory->cd( ".." ); 238 gDirectory->cd( ".." );
239 } 239 }
240 else 240 else
241 { 241 {
242 CreateHisto( histoVector, data.impl, f 242 CreateHisto( histoVector, data.impl, fullName, fullTitle, data.axes );
243 } 243 }
244 } 244 }
245 245
246 246
247 void CexmcHistoManager::CreateHisto( CexmcHis 247 void CexmcHistoManager::CreateHisto( CexmcHistoVector & histoVector,
248 CexmcHistoImpl histoI 248 CexmcHistoImpl histoImpl, const G4String & name,
249 const G4String & titl 249 const G4String & title, const CexmcHistoAxes & axes )
250 { 250 {
251 TH1 * histo( NULL ); 251 TH1 * histo( NULL );
252 252
253 switch ( histoImpl ) 253 switch ( histoImpl )
254 { 254 {
255 case Cexmc_TH1F : 255 case Cexmc_TH1F :
256 histo = new TH1F( name, title, axes.at 256 histo = new TH1F( name, title, axes.at( 0 ).nBins,
257 axes.at( 0 ).nBinsMi 257 axes.at( 0 ).nBinsMin, axes.at( 0 ).nBinsMax );
258 break; 258 break;
259 case Cexmc_TH2F : 259 case Cexmc_TH2F :
260 histo = new TH2F( name, title, axes.at 260 histo = new TH2F( name, title, axes.at( 0 ).nBins,
261 axes.at( 0 ).nBinsMi 261 axes.at( 0 ).nBinsMin, axes.at( 0 ).nBinsMax,
262 axes.at( 1 ).nBins, 262 axes.at( 1 ).nBins, axes.at( 1 ).nBinsMin,
263 axes.at( 1 ).nBinsMa 263 axes.at( 1 ).nBinsMax );
264 break; 264 break;
265 case Cexmc_TH3F : 265 case Cexmc_TH3F :
266 histo = new TH3F( name, title, axes.at 266 histo = new TH3F( name, title, axes.at( 0 ).nBins,
267 axes.at( 0 ).nBinsMi 267 axes.at( 0 ).nBinsMin, axes.at( 0 ).nBinsMax,
268 axes.at( 1 ).nBins, 268 axes.at( 1 ).nBins, axes.at( 1 ).nBinsMin,
269 axes.at( 1 ).nBinsMa 269 axes.at( 1 ).nBinsMax, axes.at( 2 ).nBins,
270 axes.at( 2 ).nBinsMi 270 axes.at( 2 ).nBinsMin, axes.at( 2 ).nBinsMax );
271 break; 271 break;
272 default : 272 default :
273 break; 273 break;
274 } 274 }
275 275
276 if ( histo ) 276 if ( histo )
277 histoVector.push_back( histo ); 277 histoVector.push_back( histo );
278 } 278 }
279 279
280 280
281 void CexmcHistoManager::Initialize( void ) 281 void CexmcHistoManager::Initialize( void )
282 { 282 {
283 if ( isInitialized ) 283 if ( isInitialized )
284 return; 284 return;
285 285
286 CexmcRunManager * runManager( static 286 CexmcRunManager * runManager( static_cast< CexmcRunManager * >(
287 G4 287 G4RunManager::GetRunManager() ) );
288 CexmcPhysicsManager * physicsManager( ru 288 CexmcPhysicsManager * physicsManager( runManager->GetPhysicsManager() );
289 289
290 if ( ! physicsManager ) 290 if ( ! physicsManager )
291 throw CexmcException( CexmcWeirdExcept 291 throw CexmcException( CexmcWeirdException );
292 292
293 CexmcProductionModel * productionModel( p 293 CexmcProductionModel * productionModel( physicsManager->
294 294 GetProductionModel() );
295 295
296 if ( ! productionModel ) 296 if ( ! productionModel )
297 throw CexmcException( CexmcWeirdExcept 297 throw CexmcException( CexmcWeirdException );
298 298
299 G4ParticleDefinition * outputParticle( 299 G4ParticleDefinition * outputParticle(
300 productionMode 300 productionModel->GetOutputParticle() );
301 G4ParticleDefinition * nucleusOutputParti 301 G4ParticleDefinition * nucleusOutputParticle(
302 productionMode 302 productionModel->GetNucleusOutputParticle() );
303 303
304 if ( ! outputParticle || ! nucleusOutputPa 304 if ( ! outputParticle || ! nucleusOutputParticle )
305 throw CexmcException( CexmcIncompleteP 305 throw CexmcException( CexmcIncompleteProductionModel );
306 306
307 opName = outputParticle->GetParticleName() 307 opName = outputParticle->GetParticleName();
308 nopName = nucleusOutputParticle->GetPartic 308 nopName = nucleusOutputParticle->GetParticleName();
309 opMass = outputParticle->GetPDGMass(); 309 opMass = outputParticle->GetPDGMass();
310 nopMass = nucleusOutputParticle->GetPDGMas 310 nopMass = nucleusOutputParticle->GetPDGMass();
311 311
312 G4String title; 312 G4String title;
313 Int_t nBinsX; 313 Int_t nBinsX;
314 Int_t nBinsY; 314 Int_t nBinsY;
315 Double_t nBinsMinX; 315 Double_t nBinsMinX;
316 Double_t nBinsMaxX; 316 Double_t nBinsMaxX;
317 Double_t nBinsMinY; 317 Double_t nBinsMinY;
318 Double_t nBinsMaxY; 318 Double_t nBinsMaxY;
319 CexmcHistoAxes axes; 319 CexmcHistoAxes axes;
320 320
321 if ( runManager->ProjectIsSaved() ) 321 if ( runManager->ProjectIsSaved() )
322 { 322 {
323 G4String projectsDir( runManager->Get 323 G4String projectsDir( runManager->GetProjectsDir() );
324 G4String resultsFile( projectsDir + " 324 G4String resultsFile( projectsDir + "/" + runManager->GetProjectId() +
325 ".root" ); 325 ".root" );
326 outFile = new TFile( resultsFile, "rec 326 outFile = new TFile( resultsFile, "recreate" );
327 } 327 }
328 else 328 else
329 { 329 {
330 outFile = new TDirectoryFile( CexmcHis 330 outFile = new TDirectoryFile( CexmcHistoDirectoryHandle,
331 CexmcHis 331 CexmcHistoDirectoryTitle );
332 gDirectory->cd( CexmcHistoDirectoryHan 332 gDirectory->cd( CexmcHistoDirectoryHandle );
333 } 333 }
334 334
335 if ( ! outFile ) 335 if ( ! outFile )
336 throw CexmcException( CexmcWeirdExcept 336 throw CexmcException( CexmcWeirdException );
337 337
338 const CexmcSetup * setup( static_cast< co 338 const CexmcSetup * setup( static_cast< const CexmcSetup * >(
339 runManager->Ge 339 runManager->GetUserDetectorConstruction() ) );
340 if ( ! setup ) 340 if ( ! setup )
341 throw CexmcException( CexmcWeirdExcept 341 throw CexmcException( CexmcWeirdException );
342 342
343 const G4LogicalVolume * lVolume( setup->G 343 const G4LogicalVolume * lVolume( setup->GetVolume( CexmcSetup::Monitor ) );
344 344
345 if ( ! lVolume ) 345 if ( ! lVolume )
346 throw CexmcException( CexmcIncompatibl 346 throw CexmcException( CexmcIncompatibleGeometry );
347 347
348 nBinsMinX = CexmcHistoBeamMomentumMin; 348 nBinsMinX = CexmcHistoBeamMomentumMin;
349 nBinsMaxX = CexmcHistoBeamMomentumMax; 349 nBinsMaxX = CexmcHistoBeamMomentumMax;
350 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 350 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) /
351 CexmcHistoBeamMomentumReso 351 CexmcHistoBeamMomentumResolution );
352 axes.push_back( CexmcHistoAxisData( nBinsX 352 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
353 AddHisto( CexmcHistoData( CexmcMomentumBP_ 353 AddHisto( CexmcHistoData( CexmcMomentumBP_TPT_Histo, Cexmc_TH1F, false,
354 false, CexmcTPT, "mombp", "Beam mome 354 false, CexmcTPT, "mombp", "Beam momentum at the monitor", axes ) );
355 AddHisto( CexmcHistoData( CexmcMomentumBP_ 355 AddHisto( CexmcHistoData( CexmcMomentumBP_RT_Histo, Cexmc_TH1F, false,
356 false, CexmcRT, "mombp", "Beam momen 356 false, CexmcRT, "mombp", "Beam momentum at the monitor", axes ) );
357 if ( verboseLevel > 0 ) 357 if ( verboseLevel > 0 )
358 { 358 {
359 AddHisto( CexmcHistoData( CexmcMomentu 359 AddHisto( CexmcHistoData( CexmcMomentumIP_TPT_Histo, Cexmc_TH1F, false,
360 false, CexmcTPT, "momip", "Momentu 360 false, CexmcTPT, "momip", "Momentum of the incident particle",
361 axes ) ); 361 axes ) );
362 } 362 }
363 363
364 G4Box * box( dynamic_cast< G4Box * >( lV 364 G4Box * box( dynamic_cast< G4Box * >( lVolume->GetSolid() ) );
365 365
366 if ( ! box ) 366 if ( ! box )
367 throw CexmcException( CexmcIncompatibl 367 throw CexmcException( CexmcIncompatibleGeometry );
368 368
369 G4double width( box->GetXHalfLength() * 2 369 G4double width( box->GetXHalfLength() * 2 );
370 G4double height( box->GetYHalfLength() * 370 G4double height( box->GetYHalfLength() * 2 );
371 G4double halfWidth( width / 2 + CexmcHist 371 G4double halfWidth( width / 2 + CexmcHistoTPSafetyArea );
372 G4double halfHeight( height / 2 + CexmcHi 372 G4double halfHeight( height / 2 + CexmcHistoTPSafetyArea );
373 373
374 nBinsX = Int_t( halfWidth * 2 / CexmcHisto 374 nBinsX = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
375 nBinsY = Int_t( halfHeight * 2 / CexmcHist 375 nBinsY = Int_t( halfHeight * 2 / CexmcHistoTPResolution );
376 axes.clear(); 376 axes.clear();
377 axes.push_back( CexmcHistoAxisData( nBinsX 377 axes.push_back( CexmcHistoAxisData( nBinsX, -halfWidth, halfWidth ) );
378 axes.push_back( CexmcHistoAxisData( nBinsY 378 axes.push_back( CexmcHistoAxisData( nBinsY, -halfHeight, halfHeight ) );
379 AddHisto( CexmcHistoData( CexmcTPInMonitor 379 AddHisto( CexmcHistoData( CexmcTPInMonitor_TPT_Histo, Cexmc_TH2F, false,
380 false, CexmcTPT, "tpmon", "Track point 380 false, CexmcTPT, "tpmon", "Track points (mon)", axes ) );
381 381
382 lVolume = setup->GetVolume( CexmcSetup::Ta 382 lVolume = setup->GetVolume( CexmcSetup::Target );
383 G4Tubs * tube( dynamic_cast< G4Tubs * >( 383 G4Tubs * tube( dynamic_cast< G4Tubs * >( lVolume->GetSolid() ) );
384 384
385 if ( ! tube ) 385 if ( ! tube )
386 throw CexmcException( CexmcIncompatibl 386 throw CexmcException( CexmcIncompatibleGeometry );
387 387
388 G4double radius( tube->GetOuterRadius() ) 388 G4double radius( tube->GetOuterRadius() );
389 height = tube->GetZHalfLength() * 2; 389 height = tube->GetZHalfLength() * 2;
390 halfWidth = radius + CexmcHistoTPSafetyAre 390 halfWidth = radius + CexmcHistoTPSafetyArea;
391 halfHeight = height / 2 + CexmcHistoTPSafe 391 halfHeight = height / 2 + CexmcHistoTPSafetyArea;
392 392
393 nBinsX = Int_t( halfWidth * 2 / CexmcHisto 393 nBinsX = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
394 nBinsY = Int_t( halfWidth * 2 / CexmcHisto 394 nBinsY = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
395 Int_t nBinsZ( Int_t( halfHeight * 2 / Cex 395 Int_t nBinsZ( Int_t( halfHeight * 2 / CexmcHistoTPResolution ) );
396 axes.clear(); 396 axes.clear();
397 axes.push_back( CexmcHistoAxisData( nBinsX 397 axes.push_back( CexmcHistoAxisData( nBinsX, -halfWidth, halfWidth ) );
398 axes.push_back( CexmcHistoAxisData( nBinsY 398 axes.push_back( CexmcHistoAxisData( nBinsY, -halfWidth, halfWidth ) );
399 axes.push_back( CexmcHistoAxisData( nBinsZ 399 axes.push_back( CexmcHistoAxisData( nBinsZ, -halfHeight, halfHeight ) );
400 AddHisto( CexmcHistoData( CexmcTPInTarget_ 400 AddHisto( CexmcHistoData( CexmcTPInTarget_TPT_Histo, Cexmc_TH3F, false,
401 false, CexmcTPT, "tptar", "Track point 401 false, CexmcTPT, "tptar", "Track points (tar)", axes ) );
402 AddHisto( CexmcHistoData( CexmcTPInTarget_ 402 AddHisto( CexmcHistoData( CexmcTPInTarget_RT_Histo, Cexmc_TH3F, false,
403 false, CexmcRT, "tptar", "Track points 403 false, CexmcRT, "tptar", "Track points (tar)", axes ) );
404 404
405 title = "Reconstructed masses (" + nopName 405 title = "Reconstructed masses (" + nopName + " vs. " + opName + ")";
406 nBinsMinX = opMass / 2; 406 nBinsMinX = opMass / 2;
407 nBinsMaxX = opMass + opMass / 2; 407 nBinsMaxX = opMass + opMass / 2;
408 nBinsMinY = nopMass / 2; 408 nBinsMinY = nopMass / 2;
409 nBinsMaxY = nopMass + nopMass / 2; 409 nBinsMaxY = nopMass + nopMass / 2;
410 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 410 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoMassResolution );
411 nBinsY = Int_t( ( nBinsMaxY - nBinsMinY ) 411 nBinsY = Int_t( ( nBinsMaxY - nBinsMinY ) / CexmcHistoMassResolution );
412 axes.clear(); 412 axes.clear();
413 axes.push_back( CexmcHistoAxisData( nBinsX 413 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
414 axes.push_back( CexmcHistoAxisData( nBinsY 414 axes.push_back( CexmcHistoAxisData( nBinsY, nBinsMinY, nBinsMaxY ) );
415 AddHisto( CexmcHistoData( CexmcRecMasses_E 415 AddHisto( CexmcHistoData( CexmcRecMasses_EDT_Histo, Cexmc_TH2F, false,
416 false, CexmcEDT, "recmasses", title, a 416 false, CexmcEDT, "recmasses", title, axes ) );
417 AddHisto( CexmcHistoData( CexmcRecMasses_R 417 AddHisto( CexmcHistoData( CexmcRecMasses_RT_Histo, Cexmc_TH2F, false,
418 false, CexmcRT, "recmasses", title, ax 418 false, CexmcRT, "recmasses", title, axes ) );
419 419
420 nBinsMinX = 0.; 420 nBinsMinX = 0.;
421 nBinsMaxX = CexmcHistoEnergyMax; 421 nBinsMaxX = CexmcHistoEnergyMax;
422 nBinsMinY = 0.; 422 nBinsMinY = 0.;
423 nBinsMaxY = CexmcHistoEnergyMax; 423 nBinsMaxY = CexmcHistoEnergyMax;
424 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 424 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoEnergyResolution );
425 nBinsY = Int_t( ( nBinsMaxY - nBinsMinY ) 425 nBinsY = Int_t( ( nBinsMaxY - nBinsMinY ) / CexmcHistoEnergyResolution );
426 axes.clear(); 426 axes.clear();
427 axes.push_back( CexmcHistoAxisData( nBinsX 427 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
428 axes.push_back( CexmcHistoAxisData( nBinsY 428 axes.push_back( CexmcHistoAxisData( nBinsY, nBinsMinY, nBinsMaxY ) );
429 AddHisto( CexmcHistoData( CexmcAbsorbedEne 429 AddHisto( CexmcHistoData( CexmcAbsorbedEnergy_EDT_Histo, Cexmc_TH2F, false,
430 false, CexmcEDT, "ae", "Absorbed energ 430 false, CexmcEDT, "ae", "Absorbed energy (rc vs. lc)", axes ) );
431 AddHisto( CexmcHistoData( CexmcAbsorbedEne 431 AddHisto( CexmcHistoData( CexmcAbsorbedEnergy_RT_Histo, Cexmc_TH2F, false,
432 false, CexmcRT, "ae", "Absorbed energy 432 false, CexmcRT, "ae", "Absorbed energy (rc vs. lc)", axes ) );
433 433
434 SetupARHistos( runManager->GetPhysicsManag 434 SetupARHistos( runManager->GetPhysicsManager()->GetProductionModel()->
435 GetAngularRanges() ); 435 GetAngularRanges() );
436 436
437 isInitialized = true; 437 isInitialized = true;
438 } 438 }
439 439
440 440
441 void CexmcHistoManager::SetupARHistos( const 441 void CexmcHistoManager::SetupARHistos( const CexmcAngularRangeList & aRanges )
442 { 442 {
443 TIter objs( gDirectory->GetList() ); 443 TIter objs( gDirectory->GetList() );
444 TObject * obj( NULL ); 444 TObject * obj( NULL );
445 445
446 while ( ( obj = ( TObject * )objs() ) ) 446 while ( ( obj = ( TObject * )objs() ) )
447 { 447 {
448 TString name( obj->GetName() ); 448 TString name( obj->GetName() );
449 449
450 if ( obj->IsFolder() && ( name.Contain 450 if ( obj->IsFolder() && ( name.Contains( TString( "_arreal_" ) ) ||
451 name.Contain 451 name.Contains( TString( "_arrec_" ) ) ) )
452 { 452 {
453 gDirectory->cd( name ); 453 gDirectory->cd( name );
454 gDirectory->DeleteAll(); 454 gDirectory->DeleteAll();
455 gDirectory->cd( ".." ); 455 gDirectory->cd( ".." );
456 } 456 }
457 } 457 }
458 458
459 for ( CexmcHistosMap::iterator k( histos. 459 for ( CexmcHistosMap::iterator k( histos.begin() ); k != histos.end();
460 460 ++k )
461 { 461 {
462 if ( k->second.empty() ) 462 if ( k->second.empty() )
463 continue; 463 continue;
464 464
465 if ( k->first >= CexmcHistoType_ARReal 465 if ( k->first >= CexmcHistoType_ARReal_START &&
466 k->first <= CexmcHistoType_ARReal 466 k->first <= CexmcHistoType_ARReal_END )
467 { 467 {
468 k->second.clear(); 468 k->second.clear();
469 } 469 }
470 } 470 }
471 471
472 for ( CexmcAngularRangeList::const_iterato 472 for ( CexmcAngularRangeList::const_iterator k( aRanges.begin() );
473 473 k != aRanges.end(); ++k )
474 { 474 {
475 AddARHistos( *k ); 475 AddARHistos( *k );
476 } 476 }
477 } 477 }
478 478
479 479
480 void CexmcHistoManager::AddARHistos( const Ce 480 void CexmcHistoManager::AddARHistos( const CexmcAngularRange & aRange )
481 { 481 {
482 G4String title; 482 G4String title;
483 Int_t nBinsX; 483 Int_t nBinsX;
484 Double_t nBinsMinX; 484 Double_t nBinsMinX;
485 Double_t nBinsMaxX; 485 Double_t nBinsMaxX;
486 CexmcHistoAxes axes; 486 CexmcHistoAxes axes;
487 487
488 title = "Reconstructed mass of " + opName; 488 title = "Reconstructed mass of " + opName;
489 nBinsMinX = opMass / 2; 489 nBinsMinX = opMass / 2;
490 nBinsMaxX = opMass + opMass / 2; 490 nBinsMaxX = opMass + opMass / 2;
491 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 491 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoMassResolution );
492 axes.push_back( CexmcHistoAxisData( nBinsX 492 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
493 AddHisto( CexmcHistoData( CexmcRecMassOP_A 493 AddHisto( CexmcHistoData( CexmcRecMassOP_ARReal_RT_Histo, Cexmc_TH1F, true,
494 false, CexmcRT, "recmassop", title, ax 494 false, CexmcRT, "recmassop", title, axes ), aRange );
495 495
496 title = "Reconstructed mass of " + nopName 496 title = "Reconstructed mass of " + nopName;
497 nBinsMinX = nopMass / 2; 497 nBinsMinX = nopMass / 2;
498 nBinsMaxX = nopMass + nopMass / 2; 498 nBinsMaxX = nopMass + nopMass / 2;
499 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 499 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoMassResolution );
500 axes.clear(); 500 axes.clear();
501 axes.push_back( CexmcHistoAxisData( nBinsX 501 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
502 AddHisto( CexmcHistoData( CexmcRecMassNOP_ 502 AddHisto( CexmcHistoData( CexmcRecMassNOP_ARReal_RT_Histo, Cexmc_TH1F, true,
503 false, CexmcRT, "recmassnop", title, a 503 false, CexmcRT, "recmassnop", title, axes ), aRange );
504 504
505 G4RunManager * runManager( G4RunManag 505 G4RunManager * runManager( G4RunManager::GetRunManager() );
506 const CexmcSetup * setup( static_cast< co 506 const CexmcSetup * setup( static_cast< const CexmcSetup * >(
507 runManager->Ge 507 runManager->GetUserDetectorConstruction() ) );
508 if ( ! setup ) 508 if ( ! setup )
509 throw CexmcException( CexmcWeirdExcept 509 throw CexmcException( CexmcWeirdException );
510 510
511 const G4LogicalVolume * lVolume( setup->G 511 const G4LogicalVolume * lVolume( setup->GetVolume(
512 512 CexmcSetup::Calorimeter ) );
513 513
514 G4Box * box( dynamic_cast< G4Box * >( lV 514 G4Box * box( dynamic_cast< G4Box * >( lVolume->GetSolid() ) );
515 515
516 if ( ! box ) 516 if ( ! box )
517 throw CexmcException( CexmcIncompatibl 517 throw CexmcException( CexmcIncompatibleGeometry );
518 518
519 G4double width( box->GetXHalfLength() * 2 519 G4double width( box->GetXHalfLength() * 2 );
520 G4double height( box->GetYHalfLength() * 520 G4double height( box->GetYHalfLength() * 2 );
521 G4double halfWidth( width / 2 + CexmcHist 521 G4double halfWidth( width / 2 + CexmcHistoTPSafetyArea );
522 G4double halfHeight( height / 2 + CexmcHi 522 G4double halfHeight( height / 2 + CexmcHistoTPSafetyArea );
523 523
524 nBinsX = Int_t( halfWidth * 2 / CexmcHisto 524 nBinsX = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
525 Int_t nBinsY( Int_t( halfHeight * 2 / Cex 525 Int_t nBinsY( Int_t( halfHeight * 2 / CexmcHistoTPResolution ) );
526 axes.clear(); 526 axes.clear();
527 axes.push_back( CexmcHistoAxisData( nBinsX 527 axes.push_back( CexmcHistoAxisData( nBinsX, -halfWidth, halfWidth ) );
528 axes.push_back( CexmcHistoAxisData( nBinsY 528 axes.push_back( CexmcHistoAxisData( nBinsY, -halfHeight, halfHeight ) );
529 529
530 /* looks like there is no possibility to d 530 /* looks like there is no possibility to draw descending xaxis in ROOT,
531 * so imagine that you look at calorimeter 531 * so imagine that you look at calorimeters from behind, i.e. your face to
532 * the beam */ 532 * the beam */
533 AddHisto( CexmcHistoData( CexmcOPDPAtLeftC 533 AddHisto( CexmcHistoData( CexmcOPDPAtLeftCalorimeter_ARReal_EDT_Histo,
534 Cexmc_TH2F, true, false, CexmcEDT, "op 534 Cexmc_TH2F, true, false, CexmcEDT, "opdpcl",
535 "Gamma position on the surface (lc)", 535 "Gamma position on the surface (lc)", axes ), aRange );
536 AddHisto( CexmcHistoData( CexmcOPDPAtRight 536 AddHisto( CexmcHistoData( CexmcOPDPAtRightCalorimeter_ARReal_EDT_Histo,
537 Cexmc_TH2F, true, false, CexmcEDT, "op 537 Cexmc_TH2F, true, false, CexmcEDT, "opdpcr",
538 "Gamma position on the surface (rc)", 538 "Gamma position on the surface (rc)", axes ), aRange );
539 AddHisto( CexmcHistoData( CexmcOPDPAtLeftC 539 AddHisto( CexmcHistoData( CexmcOPDPAtLeftCalorimeter_ARReal_RT_Histo,
540 Cexmc_TH2F, true, false, CexmcRT, "opd 540 Cexmc_TH2F, true, false, CexmcRT, "opdpcl",
541 "Gamma position on the surface (lc)", 541 "Gamma position on the surface (lc)", axes ), aRange );
542 AddHisto( CexmcHistoData( CexmcOPDPAtRight 542 AddHisto( CexmcHistoData( CexmcOPDPAtRightCalorimeter_ARReal_RT_Histo,
543 Cexmc_TH2F, true, false, CexmcRT, "opd 543 Cexmc_TH2F, true, false, CexmcRT, "opdpcr",
544 "Gamma position on the surface (rc)", 544 "Gamma position on the surface (rc)", axes ), aRange );
545 AddHisto( CexmcHistoData( CexmcRecOPDPAtLe 545 AddHisto( CexmcHistoData( CexmcRecOPDPAtLeftCalorimeter_ARReal_EDT_Histo,
546 Cexmc_TH2F, true, false, CexmcEDT, "re 546 Cexmc_TH2F, true, false, CexmcEDT, "recopdpcl",
547 "Reconstructed gamma position on the s 547 "Reconstructed gamma position on the surface (lc)", axes ), aRange );
548 AddHisto( CexmcHistoData( CexmcRecOPDPAtRi 548 AddHisto( CexmcHistoData( CexmcRecOPDPAtRightCalorimeter_ARReal_EDT_Histo,
549 Cexmc_TH2F, true, false, CexmcEDT, "re 549 Cexmc_TH2F, true, false, CexmcEDT, "recopdpcr",
550 "Reconstructed gamma position on the s 550 "Reconstructed gamma position on the surface (rc)", axes ), aRange );
551 AddHisto( CexmcHistoData( CexmcRecOPDPAtLe 551 AddHisto( CexmcHistoData( CexmcRecOPDPAtLeftCalorimeter_ARReal_RT_Histo,
552 Cexmc_TH2F, true, false, CexmcRT, "rec 552 Cexmc_TH2F, true, false, CexmcRT, "recopdpcl",
553 "Reconstructed gamma position on the s 553 "Reconstructed gamma position on the surface (lc)", axes ), aRange );
554 AddHisto( CexmcHistoData( CexmcRecOPDPAtRi 554 AddHisto( CexmcHistoData( CexmcRecOPDPAtRightCalorimeter_ARReal_RT_Histo,
555 Cexmc_TH2F, true, false, CexmcRT, "rec 555 Cexmc_TH2F, true, false, CexmcRT, "recopdpcr",
556 "Reconstructed gamma position on the s 556 "Reconstructed gamma position on the surface (rc)", axes ), aRange );
557 557
558 nBinsMinX = 0.; 558 nBinsMinX = 0.;
559 nBinsMaxX = CexmcHistoEnergyMax; 559 nBinsMaxX = CexmcHistoEnergyMax;
560 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 560 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoEnergyResolution );
561 axes.clear(); 561 axes.clear();
562 axes.push_back( CexmcHistoAxisData( nBinsX 562 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
563 AddHisto( CexmcHistoData( CexmcKinEnAtLeft 563 AddHisto( CexmcHistoData( CexmcKinEnAtLeftCalorimeter_ARReal_TPT_Histo,
564 Cexmc_TH1F, true, false, CexmcTPT, "ke 564 Cexmc_TH1F, true, false, CexmcTPT, "kecl",
565 "Kinetic energy of gamma (lc)", axes ) 565 "Kinetic energy of gamma (lc)", axes ), aRange );
566 AddHisto( CexmcHistoData( CexmcKinEnAtRigh 566 AddHisto( CexmcHistoData( CexmcKinEnAtRightCalorimeter_ARReal_TPT_Histo,
567 Cexmc_TH1F, true, false, CexmcTPT, "ke 567 Cexmc_TH1F, true, false, CexmcTPT, "kecr",
568 "Kinetic energy of gamma (rc)", axes ) 568 "Kinetic energy of gamma (rc)", axes ), aRange );
569 AddHisto( CexmcHistoData( CexmcKinEnAtLeft 569 AddHisto( CexmcHistoData( CexmcKinEnAtLeftCalorimeter_ARReal_RT_Histo,
570 Cexmc_TH1F, true, false, CexmcRT, "kec 570 Cexmc_TH1F, true, false, CexmcRT, "kecl",
571 "Kinetic energy of gamma (lc)", axes ) 571 "Kinetic energy of gamma (lc)", axes ), aRange );
572 AddHisto( CexmcHistoData( CexmcKinEnAtRigh 572 AddHisto( CexmcHistoData( CexmcKinEnAtRightCalorimeter_ARReal_RT_Histo,
573 Cexmc_TH1F, true, false, CexmcRT, "kec 573 Cexmc_TH1F, true, false, CexmcRT, "kecr",
574 "Kinetic energy of gamma (rc)", axes ) 574 "Kinetic energy of gamma (rc)", axes ), aRange );
575 AddHisto( CexmcHistoData( CexmcAbsEnInLeft 575 AddHisto( CexmcHistoData( CexmcAbsEnInLeftCalorimeter_ARReal_EDT_Histo,
576 Cexmc_TH1F, true, false, CexmcEDT, "ae 576 Cexmc_TH1F, true, false, CexmcEDT, "aecl",
577 "Absorbed energy (lc)", axes ), aRange 577 "Absorbed energy (lc)", axes ), aRange );
578 AddHisto( CexmcHistoData( CexmcAbsEnInRigh 578 AddHisto( CexmcHistoData( CexmcAbsEnInRightCalorimeter_ARReal_EDT_Histo,
579 Cexmc_TH1F, true, false, CexmcEDT, "ae 579 Cexmc_TH1F, true, false, CexmcEDT, "aecr",
580 "Absorbed energy (rc)", axes ), aRange 580 "Absorbed energy (rc)", axes ), aRange );
581 AddHisto( CexmcHistoData( CexmcAbsEnInLeft 581 AddHisto( CexmcHistoData( CexmcAbsEnInLeftCalorimeter_ARReal_RT_Histo,
582 Cexmc_TH1F, true, false, CexmcRT, "aec 582 Cexmc_TH1F, true, false, CexmcRT, "aecl",
583 "Absorbed energy (lc)", axes ), aRange 583 "Absorbed energy (lc)", axes ), aRange );
584 AddHisto( CexmcHistoData( CexmcAbsEnInRigh 584 AddHisto( CexmcHistoData( CexmcAbsEnInRightCalorimeter_ARReal_RT_Histo,
585 Cexmc_TH1F, true, false, CexmcRT, "aec 585 Cexmc_TH1F, true, false, CexmcRT, "aecr",
586 "Absorbed energy (rc)", axes ), aRange 586 "Absorbed energy (rc)", axes ), aRange );
587 587
588 nBinsMinX = CexmcHistoMissEnergyMin; 588 nBinsMinX = CexmcHistoMissEnergyMin;
589 nBinsMaxX = CexmcHistoMissEnergyMax; 589 nBinsMaxX = CexmcHistoMissEnergyMax;
590 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 590 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) /
591 CexmcHistoMissEnergyResolu 591 CexmcHistoMissEnergyResolution );
592 axes.clear(); 592 axes.clear();
593 axes.push_back( CexmcHistoAxisData( nBinsX 593 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
594 AddHisto( CexmcHistoData( CexmcMissEnFromL 594 AddHisto( CexmcHistoData( CexmcMissEnFromLeftCalorimeter_ARReal_RT_Histo,
595 Cexmc_TH1F, true, false, CexmcRT, "mec 595 Cexmc_TH1F, true, false, CexmcRT, "mecl",
596 "Missing energy (lc)", axes ), aRange 596 "Missing energy (lc)", axes ), aRange );
597 AddHisto( CexmcHistoData( CexmcMissEnFromR 597 AddHisto( CexmcHistoData( CexmcMissEnFromRightCalorimeter_ARReal_RT_Histo,
598 Cexmc_TH1F, true, false, CexmcRT, "mec 598 Cexmc_TH1F, true, false, CexmcRT, "mecr",
599 "Missing energy (rc)", axes ), aRange 599 "Missing energy (rc)", axes ), aRange );
600 600
601 title = "Kinetic energy of newborn " + opN 601 title = "Kinetic energy of newborn " + opName + " (lab)";
602 nBinsMinX = 0.; 602 nBinsMinX = 0.;
603 nBinsMaxX = CexmcHistoEnergyMax; 603 nBinsMaxX = CexmcHistoEnergyMax;
604 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 604 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoEnergyResolution );
605 axes.clear(); 605 axes.clear();
606 axes.push_back( CexmcHistoAxisData( nBinsX 606 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
607 AddHisto( CexmcHistoData( CexmcKinEnOP_LAB 607 AddHisto( CexmcHistoData( CexmcKinEnOP_LAB_ARReal_TPT_Histo,
608 Cexmc_TH1F, true, false, CexmcTPT, "ke 608 Cexmc_TH1F, true, false, CexmcTPT, "keop_lab", title, axes ), aRange );
609 AddHisto( CexmcHistoData( CexmcKinEnOP_LAB 609 AddHisto( CexmcHistoData( CexmcKinEnOP_LAB_ARReal_RT_Histo,
610 Cexmc_TH1F, true, false, CexmcRT, "keo 610 Cexmc_TH1F, true, false, CexmcRT, "keop_lab", title, axes ), aRange );
611 611
612 title = "Angle of newborn " + opName + " ( 612 title = "Angle of newborn " + opName + " (scm)";
613 nBinsMinX = -1.0; 613 nBinsMinX = -1.0;
614 nBinsMaxX = 1.0; 614 nBinsMaxX = 1.0;
615 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 615 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoAngularCResolution );
616 axes.clear(); 616 axes.clear();
617 axes.push_back( CexmcHistoAxisData( nBinsX 617 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
618 AddHisto( CexmcHistoData( CexmcAngleOP_SCM 618 AddHisto( CexmcHistoData( CexmcAngleOP_SCM_ARReal_TPT_Histo,
619 Cexmc_TH1F, true, false, CexmcTPT, "ao 619 Cexmc_TH1F, true, false, CexmcTPT, "aop_scm", title, axes ), aRange );
620 AddHisto( CexmcHistoData( CexmcAngleOP_SCM 620 AddHisto( CexmcHistoData( CexmcAngleOP_SCM_ARReal_RT_Histo,
621 Cexmc_TH1F, true, false, CexmcRT, "aop 621 Cexmc_TH1F, true, false, CexmcRT, "aop_scm", title, axes ), aRange );
622 622
623 title = "Reconstruced angle of newborn " + 623 title = "Reconstruced angle of newborn " + opName + " (scm)";
624 AddHisto( CexmcHistoData( CexmcRecAngleOP_ 624 AddHisto( CexmcHistoData( CexmcRecAngleOP_SCM_ARReal_RT_Histo,
625 Cexmc_TH1F, true, false, CexmcRT, "rec 625 Cexmc_TH1F, true, false, CexmcRT, "recaop_scm", title, axes ), aRange );
626 626
627 title = "Real - reconstruced angle of newb 627 title = "Real - reconstruced angle of newborn " + opName + " (scm)";
628 AddHisto( CexmcHistoData( CexmcDiffAngleOP 628 AddHisto( CexmcHistoData( CexmcDiffAngleOP_SCM_ARReal_RT_Histo,
629 Cexmc_TH1F, true, false, CexmcRT, "dif 629 Cexmc_TH1F, true, false, CexmcRT, "diffaop_scm", title, axes ),
630 aRange ); 630 aRange );
631 631
632 nBinsMinX = 0.; 632 nBinsMinX = 0.;
633 nBinsMaxX = 360.; 633 nBinsMaxX = 360.;
634 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 634 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoAngularResolution );
635 axes.clear(); 635 axes.clear();
636 axes.push_back( CexmcHistoAxisData( nBinsX 636 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
637 AddHisto( CexmcHistoData( CexmcOpenAngle_A 637 AddHisto( CexmcHistoData( CexmcOpenAngle_ARReal_TPT_Histo,
638 Cexmc_TH1F, true, false, CexmcTPT, "oa 638 Cexmc_TH1F, true, false, CexmcTPT, "oa",
639 "Open angle between the gammas", axes 639 "Open angle between the gammas", axes ), aRange );
640 AddHisto( CexmcHistoData( CexmcOpenAngle_A 640 AddHisto( CexmcHistoData( CexmcOpenAngle_ARReal_RT_Histo,
641 Cexmc_TH1F, true, false, CexmcRT, "oa" 641 Cexmc_TH1F, true, false, CexmcRT, "oa",
642 "Open angle between the gammas", axes 642 "Open angle between the gammas", axes ), aRange );
643 AddHisto( CexmcHistoData( CexmcRecOpenAngl 643 AddHisto( CexmcHistoData( CexmcRecOpenAngle_ARReal_RT_Histo,
644 Cexmc_TH1F, true, false, CexmcRT, "rec 644 Cexmc_TH1F, true, false, CexmcRT, "recoa",
645 "Reconstructed open angle between the 645 "Reconstructed open angle between the gammas", axes ), aRange );
646 646
647 nBinsMinX = -180.; 647 nBinsMinX = -180.;
648 nBinsMaxX = 180.; 648 nBinsMaxX = 180.;
649 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 649 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoAngularResolution );
650 axes.clear(); 650 axes.clear();
651 axes.push_back( CexmcHistoAxisData( nBinsX 651 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
652 AddHisto( CexmcHistoData( CexmcDiffOpenAng 652 AddHisto( CexmcHistoData( CexmcDiffOpenAngle_ARReal_RT_Histo,
653 Cexmc_TH1F, true, false, CexmcRT, "dif 653 Cexmc_TH1F, true, false, CexmcRT, "diffoa",
654 "Real - reconstructed open angle betwe 654 "Real - reconstructed open angle between the gammas", axes ), aRange );
655 655
656 lVolume = setup->GetVolume( CexmcSetup::Ta 656 lVolume = setup->GetVolume( CexmcSetup::Target );
657 G4Tubs * tube( dynamic_cast< G4Tubs * >( 657 G4Tubs * tube( dynamic_cast< G4Tubs * >( lVolume->GetSolid() ) );
658 658
659 if ( ! tube ) 659 if ( ! tube )
660 throw CexmcException( CexmcIncompatibl 660 throw CexmcException( CexmcIncompatibleGeometry );
661 661
662 G4double radius( tube->GetOuterRadius() ) 662 G4double radius( tube->GetOuterRadius() );
663 height = tube->GetZHalfLength() * 2; 663 height = tube->GetZHalfLength() * 2;
664 halfWidth = radius + CexmcHistoTPSafetyAre 664 halfWidth = radius + CexmcHistoTPSafetyArea;
665 halfHeight = height / 2 + CexmcHistoTPSafe 665 halfHeight = height / 2 + CexmcHistoTPSafetyArea;
666 666
667 nBinsX = Int_t( halfWidth * 2 / CexmcHisto 667 nBinsX = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
668 nBinsY = Int_t( halfWidth * 2 / CexmcHisto 668 nBinsY = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
669 Int_t nBinsZ( Int_t( halfHeight * 2 / Cex 669 Int_t nBinsZ( Int_t( halfHeight * 2 / CexmcHistoTPResolution ) );
670 axes.clear(); 670 axes.clear();
671 axes.push_back( CexmcHistoAxisData( nBinsX 671 axes.push_back( CexmcHistoAxisData( nBinsX, -halfWidth, halfWidth ) );
672 axes.push_back( CexmcHistoAxisData( nBinsY 672 axes.push_back( CexmcHistoAxisData( nBinsY, -halfWidth, halfWidth ) );
673 axes.push_back( CexmcHistoAxisData( nBinsZ 673 axes.push_back( CexmcHistoAxisData( nBinsZ, -halfHeight, halfHeight ) );
674 AddHisto( CexmcHistoData( CexmcTPInTarget_ 674 AddHisto( CexmcHistoData( CexmcTPInTarget_ARReal_TPT_Histo, Cexmc_TH3F,
675 true, false, CexmcTPT, "tptar", "Track 675 true, false, CexmcTPT, "tptar", "Track points (tar)", axes ), aRange );
676 AddHisto( CexmcHistoData( CexmcTPInTarget_ 676 AddHisto( CexmcHistoData( CexmcTPInTarget_ARReal_RT_Histo, Cexmc_TH3F,
677 true, false, CexmcRT, "tptar", "Track 677 true, false, CexmcRT, "tptar", "Track points (tar)", axes ), aRange );
678 } 678 }
679 679
680 680
681 void CexmcHistoManager::Add( CexmcHistoType 681 void CexmcHistoManager::Add( CexmcHistoType histoType, unsigned int index,
682 G4double x ) 682 G4double x )
683 { 683 {
684 CexmcHistosMap::iterator found( histos.fi 684 CexmcHistosMap::iterator found( histos.find( histoType ) );
685 if ( found == histos.end() || histos[ hist 685 if ( found == histos.end() || histos[ histoType ].size() <= index )
686 throw CexmcException( CexmcWeirdExcept 686 throw CexmcException( CexmcWeirdException );
687 687
688 histos[ histoType ][ index ]->Fill( x ); 688 histos[ histoType ][ index ]->Fill( x );
689 } 689 }
690 690
691 691
692 void CexmcHistoManager::Add( CexmcHistoType 692 void CexmcHistoManager::Add( CexmcHistoType histoType, unsigned int index,
693 G4double x, G4d 693 G4double x, G4double y )
694 { 694 {
695 CexmcHistosMap::iterator found( histos.fi 695 CexmcHistosMap::iterator found( histos.find( histoType ) );
696 if ( found == histos.end() || histos[ hist 696 if ( found == histos.end() || histos[ histoType ].size() <= index )
697 throw CexmcException( CexmcWeirdExcept 697 throw CexmcException( CexmcWeirdException );
698 698
699 /* no cast needed because TH1 has virtual 699 /* no cast needed because TH1 has virtual method
700 * Fill( Double_t, Double_t ) */ 700 * Fill( Double_t, Double_t ) */
701 histos[ histoType ][ index ]->Fill( x, y ) 701 histos[ histoType ][ index ]->Fill( x, y );
702 } 702 }
703 703
704 704
705 void CexmcHistoManager::Add( CexmcHistoType 705 void CexmcHistoManager::Add( CexmcHistoType histoType, unsigned int index,
706 G4double x, G4d 706 G4double x, G4double y, G4double z )
707 { 707 {
708 CexmcHistosMap::iterator found( histos.fi 708 CexmcHistosMap::iterator found( histos.find( histoType ) );
709 if ( found == histos.end() || histos[ hist 709 if ( found == histos.end() || histos[ histoType ].size() <= index )
710 throw CexmcException( CexmcWeirdExcept 710 throw CexmcException( CexmcWeirdException );
711 711
712 /* cast needed because TH1 does not have v 712 /* cast needed because TH1 does not have virtual method
713 * Fill( Double_t, Double_t, Double_t ) */ 713 * Fill( Double_t, Double_t, Double_t ) */
714 TH3 * histo( static_cast< TH3 * >( histos 714 TH3 * histo( static_cast< TH3 * >( histos[ histoType ][ index ] ) );
715 715
716 histo->Fill( x, y, z ); 716 histo->Fill( x, y, z );
717 } 717 }
718 718
719 719
720 void CexmcHistoManager::Add( CexmcHistoType 720 void CexmcHistoManager::Add( CexmcHistoType histoType, unsigned int index,
721 G4int binX, G4i 721 G4int binX, G4int binY, G4double value )
722 { 722 {
723 CexmcHistosMap::iterator found( histos.fi 723 CexmcHistosMap::iterator found( histos.find( histoType ) );
724 if ( found == histos.end() || histos[ hist 724 if ( found == histos.end() || histos[ histoType ].size() <= index )
725 throw CexmcException( CexmcWeirdExcept 725 throw CexmcException( CexmcWeirdException );
726 726
727 ++binX; 727 ++binX;
728 ++binY; 728 ++binY;
729 Double_t curValue( histos[ histoType ][ i 729 Double_t curValue( histos[ histoType ][ index ]->GetBinContent(
730 730 binX, binY ) );
731 histos[ histoType ][ index ]->SetBinConten 731 histos[ histoType ][ index ]->SetBinContent( binX, binY,
732 732 curValue + value / GeV );
733 } 733 }
734 734
735 735
736 void CexmcHistoManager::List( void ) const 736 void CexmcHistoManager::List( void ) const
737 { 737 {
738 /* BEWARE: list will be printed on stdout 738 /* BEWARE: list will be printed on stdout */
739 gDirectory->ls(); 739 gDirectory->ls();
740 } 740 }
741 741
742 742
743 void CexmcHistoManager::Print( const G4String 743 void CexmcHistoManager::Print( const G4String & value )
744 { 744 {
745 TObject * histo( gDirectory->FindObjectAn 745 TObject * histo( gDirectory->FindObjectAny( value.c_str() ) );
746 746
747 if ( ! histo ) 747 if ( ! histo )
748 { 748 {
749 G4cout << "Histogram '" << value << "' 749 G4cout << "Histogram '" << value << "' was not found" << G4endl;
750 return; 750 return;
751 } 751 }
752 752
753 /* BEWARE: histo will be printed on stdout 753 /* BEWARE: histo will be printed on stdout */
754 histo->Print( "range" ); 754 histo->Print( "range" );
755 } 755 }
756 756
757 757
758 #ifdef CEXMC_USE_ROOTQT 758 #ifdef CEXMC_USE_ROOTQT
759 759
760 void CexmcHistoManager::Draw( const G4String 760 void CexmcHistoManager::Draw( const G4String & histoName,
761 const G4String 761 const G4String & histoDrawOptions )
762 { 762 {
763 if ( ! areLiveHistogramsEnabled ) 763 if ( ! areLiveHistogramsEnabled )
764 { 764 {
765 G4cout << "Live histograms option is d 765 G4cout << "Live histograms option is disabled" << G4endl;
766 return; 766 return;
767 } 767 }
768 768
769 TObject * histo( gDirectory->FindObjectAn 769 TObject * histo( gDirectory->FindObjectAny( histoName ) );
770 770
771 if ( ! histo ) 771 if ( ! histo )
772 { 772 {
773 G4cout << "Histogram '" << histoName < 773 G4cout << "Histogram '" << histoName << "' was not found" << G4endl;
774 return; 774 return;
775 } 775 }
776 776
777 if ( ! rootCanvas ) 777 if ( ! rootCanvas )
778 { 778 {
779 /* save default application font becau 779 /* save default application font because rootCanvas will break it */
780 QFont defaultAppFont( QApplication::f 780 QFont defaultAppFont( QApplication::font() );
781 rootCanvas = new CexmcHistoWidget; 781 rootCanvas = new CexmcHistoWidget;
782 QApplication::setFont( defaultAppFont 782 QApplication::setFont( defaultAppFont );
783 rootCanvas->resize( CexmcHistoCanvasWi 783 rootCanvas->resize( CexmcHistoCanvasWidth, CexmcHistoCanvasHeight );
784 rootCanvas->GetCanvas()->cd(); 784 rootCanvas->GetCanvas()->cd();
785 } 785 }
786 786
787 histo->Draw( histoDrawOptions ); 787 histo->Draw( histoDrawOptions );
788 rootCanvas->show(); 788 rootCanvas->show();
789 rootCanvas->GetCanvas()->Update(); 789 rootCanvas->GetCanvas()->Update();
790 } 790 }
791 791
792 792
793 void CexmcHistoManager::EnableLiveHistograms( 793 void CexmcHistoManager::EnableLiveHistograms( G4UIsession * session,
794 794 G4bool on )
795 { 795 {
796 areLiveHistogramsEnabled = on; 796 areLiveHistogramsEnabled = on;
797 797
798 if ( ! on || ! isInitialized ) 798 if ( ! on || ! isInitialized )
799 return; 799 return;
800 800
801 G4UIQt * qtSession( dynamic_cast< G4UIQt 801 G4UIQt * qtSession( dynamic_cast< G4UIQt * >( session ) );
802 802
803 if ( ! qtSession ) 803 if ( ! qtSession )
804 return; 804 return;
805 805
806 if ( ! histoMenuHandle.empty() && ! isHist 806 if ( ! histoMenuHandle.empty() && ! isHistoMenuInitialized )
807 { 807 {
808 qtSession->AddMenu( histoMenuHandle, h 808 qtSession->AddMenu( histoMenuHandle, histoMenuLabel );
809 BuildMenuTree( qtSession, histoMenuHan 809 BuildMenuTree( qtSession, histoMenuHandle, gDirectory->GetList() );
810 isHistoMenuInitialized = true; 810 isHistoMenuInitialized = true;
811 } 811 }
812 } 812 }
813 813
814 814
815 void CexmcHistoManager::BuildMenuTree( G4UIQt 815 void CexmcHistoManager::BuildMenuTree( G4UIQt * session,
816 const 816 const G4String & menu, TList * ls )
817 { 817 {
818 TIter objs( ls ); 818 TIter objs( ls );
819 TObject * obj( NULL ); 819 TObject * obj( NULL );
820 820
821 while ( ( obj = ( TObject * )objs() ) ) 821 while ( ( obj = ( TObject * )objs() ) )
822 { 822 {
823 G4String name( obj->GetName() ); 823 G4String name( obj->GetName() );
824 G4String title( obj->GetTitle() ); 824 G4String title( obj->GetTitle() );
825 825
826 if ( obj->IsFolder() ) 826 if ( obj->IsFolder() )
827 { 827 {
828 AddSubmenu( session, menu, name, t 828 AddSubmenu( session, menu, name, title );
829 BuildMenuTree( session, name, ( ( 829 BuildMenuTree( session, name, ( ( TDirectory * )obj )->GetList() );
830 } 830 }
831 else 831 else
832 { 832 {
833 G4String options( name ); 833 G4String options( name );
834 834
835 do 835 do
836 { 836 {
837 if ( obj->InheritsFrom( TH3::C 837 if ( obj->InheritsFrom( TH3::Class() ) &&
838 ! drawOptions3D.empty() ) 838 ! drawOptions3D.empty() )
839 { 839 {
840 title = G4String( "3: " ) 840 title = G4String( "3: " ) + title;
841 options += G4String( " " ) 841 options += G4String( " " ) + drawOptions3D;
842 break; 842 break;
843 } 843 }
844 if ( obj->InheritsFrom( TH2::C 844 if ( obj->InheritsFrom( TH2::Class() ) &&
845 ! drawOptions2D.empty() ) 845 ! drawOptions2D.empty() )
846 { 846 {
847 title = G4String( "2: " ) 847 title = G4String( "2: " ) + title;
848 options += G4String( " " ) 848 options += G4String( " " ) + drawOptions2D;
849 break; 849 break;
850 } 850 }
851 if ( obj->InheritsFrom( TH1::C 851 if ( obj->InheritsFrom( TH1::Class() ) &&
852 ! drawOptions1D.empty() ) 852 ! drawOptions1D.empty() )
853 { 853 {
854 options += G4String( " " ) 854 options += G4String( " " ) + drawOptions1D;
855 break; 855 break;
856 } 856 }
857 } while ( false ); 857 } while ( false );
858 858
859 G4String cmd( CexmcMessenger::his 859 G4String cmd( CexmcMessenger::histoDirName + "draw " + options );
860 session->AddButton( menu, title.c_ 860 session->AddButton( menu, title.c_str(), cmd );
861 } 861 }
862 } 862 }
863 } 863 }
864 864
865 865
866 void CexmcHistoManager::AddSubmenu( G4UIQt * 866 void CexmcHistoManager::AddSubmenu( G4UIQt * session,
867 const G4S 867 const G4String & parent,
868 const G4S 868 const G4String & name,
869 const G4S 869 const G4String & label )
870 { 870 {
871 QMenu * menu( new QMenu( label.c_str() ) ); 871 QMenu * menu( new QMenu( label.c_str() ) );
872 QMenu * parentMenu( ( QMenu * )session->Get 872 QMenu * parentMenu( ( QMenu * )session->GetInteractor( parent ) );
873 873
874 parentMenu->addMenu( menu ); 874 parentMenu->addMenu( menu );
875 session->AddInteractor( name, ( G4Interactor 875 session->AddInteractor( name, ( G4Interactor )menu );
876 } 876 }
877 877
878 #endif 878 #endif
879 879
880 #endif 880 #endif
881 881
882 882