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