| 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 #include <TRegexp.h>
57 #ifdef CEXMC_USE_ROOTQT 58 #ifdef CEXMC_USE_ROOTQT
58 #include <TCanvas.h> 59 #include <TCanvas.h>
59 #include <TList.h> <<
60 #include <QApplication> 60 #include <QApplication>
61 #include <QFont> 61 #include <QFont>
62 #include <QMenu> <<
63 #include <G4UIQt.hh> <<
64 #include "CexmcMessenger.hh" <<
65 #endif 62 #endif
66 #include <G4LogicalVolume.hh> 63 #include <G4LogicalVolume.hh>
67 #include <G4Box.hh> 64 #include <G4Box.hh>
68 #include <G4Tubs.hh> 65 #include <G4Tubs.hh>
69 #include <G4SystemOfUnits.hh> 66 #include <G4SystemOfUnits.hh>
70 #include "CexmcHistoManager.hh" 67 #include "CexmcHistoManager.hh"
71 #include "CexmcHistoManagerMessenger.hh" 68 #include "CexmcHistoManagerMessenger.hh"
72 #include "CexmcProductionModel.hh" 69 #include "CexmcProductionModel.hh"
73 #include "CexmcPhysicsManager.hh" 70 #include "CexmcPhysicsManager.hh"
74 #include "CexmcRunManager.hh" 71 #include "CexmcRunManager.hh"
75 #include "CexmcSetup.hh" 72 #include "CexmcSetup.hh"
76 #include "CexmcException.hh" 73 #include "CexmcException.hh"
77 #include "CexmcHistoWidget.hh" 74 #include "CexmcHistoWidget.hh"
78 75
79 76
80 namespace 77 namespace
81 { 78 {
82 const G4double CexmcHistoBeamMomentumMin( 79 const G4double CexmcHistoBeamMomentumMin( 0.0 * GeV );
83 const G4double CexmcHistoBeamMomentumMax( 80 const G4double CexmcHistoBeamMomentumMax( 1.0 * GeV );
84 const G4double CexmcHistoBeamMomentumReso 81 const G4double CexmcHistoBeamMomentumResolution( 0.5 * MeV );
85 const G4double CexmcHistoTPResolution( 0. 82 const G4double CexmcHistoTPResolution( 0.1 * cm );
86 const G4double CexmcHistoTPSafetyArea( 1. 83 const G4double CexmcHistoTPSafetyArea( 1.0 * cm );
87 const G4double CexmcHistoMassResolution( 84 const G4double CexmcHistoMassResolution( 1.0 * MeV );
88 const G4double CexmcHistoEnergyMax( 1.0 * 85 const G4double CexmcHistoEnergyMax( 1.0 * GeV );
89 const G4double CexmcHistoEnergyResolution 86 const G4double CexmcHistoEnergyResolution( 1.0 * MeV );
90 const G4double CexmcHistoMissEnergyMin( - 87 const G4double CexmcHistoMissEnergyMin( -0.1 * GeV );
91 const G4double CexmcHistoMissEnergyMax( 0 88 const G4double CexmcHistoMissEnergyMax( 0.2 * GeV );
92 const G4double CexmcHistoMissEnergyResolu 89 const G4double CexmcHistoMissEnergyResolution( 0.2 * MeV );
93 const G4double CexmcHistoAngularResolutio 90 const G4double CexmcHistoAngularResolution( 0.5 );
94 const G4double CexmcHistoAngularCResoluti 91 const G4double CexmcHistoAngularCResolution( 0.001 );
95 #ifdef CEXMC_USE_ROOTQT <<
96 const G4int CexmcHistoCanvasWidth( 800 92 const G4int CexmcHistoCanvasWidth( 800 );
97 const G4int CexmcHistoCanvasHeight( 60 93 const G4int CexmcHistoCanvasHeight( 600 );
98 #endif <<
99 const G4String CexmcHistoDirectoryHandle( <<
100 const G4String CexmcHistoDirectoryTitle( <<
101 } 94 }
102 95
103 96
104 CexmcHistoManager * CexmcHistoManager::instan 97 CexmcHistoManager * CexmcHistoManager::instance( NULL );
105 98
106 99
107 CexmcHistoManager * CexmcHistoManager::Instan 100 CexmcHistoManager * CexmcHistoManager::Instance( void )
108 { 101 {
109 if ( instance == NULL ) 102 if ( instance == NULL )
110 instance = new CexmcHistoManager; 103 instance = new CexmcHistoManager;
111 104
112 return instance; 105 return instance;
113 } 106 }
114 107
115 108
116 void CexmcHistoManager::Destroy( void ) 109 void CexmcHistoManager::Destroy( void )
117 { 110 {
118 delete instance; 111 delete instance;
119 instance = NULL; 112 instance = NULL;
120 } 113 }
121 114
122 115
123 CexmcHistoManager::CexmcHistoManager() : outFi 116 CexmcHistoManager::CexmcHistoManager() : outFile( NULL ),
124 isInitialized( false ), opName( "" ), nopN 117 isInitialized( false ), opName( "" ), nopName( "" ), opMass( 0. ),
125 nopMass( 0. ), verboseLevel( 0 ), 118 nopMass( 0. ), verboseLevel( 0 ),
126 #ifdef CEXMC_USE_ROOTQT 119 #ifdef CEXMC_USE_ROOTQT
127 rootCanvas( NULL ), areLiveHistogramsEnabl << 120 rootCanvas( NULL ),
128 isHistoMenuInitialized( false ), drawOptio <<
129 drawOptions3D( "" ), histoMenuHandle( "" ) <<
130 #endif 121 #endif
131 messenger( NULL ) 122 messenger( NULL )
132 { 123 {
133 for ( int i( 0 ); i < CexmcHistoType_SIZE 124 for ( int i( 0 ); i < CexmcHistoType_SIZE; ++i )
134 { 125 {
135 histos.insert( CexmcHistoPair( CexmcHi 126 histos.insert( CexmcHistoPair( CexmcHistoType( i ),
136 CexmcHi 127 CexmcHistoVector() ) );
137 } 128 }
138 129
139 messenger = new CexmcHistoManagerMessenger 130 messenger = new CexmcHistoManagerMessenger( this );
140 } 131 }
141 132
142 133
143 CexmcHistoManager::~CexmcHistoManager() 134 CexmcHistoManager::~CexmcHistoManager()
144 { 135 {
145 if ( outFile ) 136 if ( outFile )
146 { 137 {
147 outFile->Write(); 138 outFile->Write();
148 outFile->Close(); 139 outFile->Close();
149 } 140 }
150 141
151 /* all histograms will be deleted by outFi 142 /* all histograms will be deleted by outFile destructor! */
152 delete outFile; 143 delete outFile;
153 #ifdef CEXMC_USE_ROOTQT 144 #ifdef CEXMC_USE_ROOTQT
154 delete rootCanvas; 145 delete rootCanvas;
155 #endif 146 #endif
156 delete messenger; 147 delete messenger;
157 } 148 }
158 149
159 150
160 void CexmcHistoManager::AddHisto( const Cexmc 151 void CexmcHistoManager::AddHisto( const CexmcHistoData & data,
161 const Cexmc 152 const CexmcAngularRange & aRange )
162 { 153 {
163 G4String fullName( data.name ); 154 G4String fullName( data.name );
164 G4String fullTitle( data.title ); 155 G4String fullTitle( data.title );
165 G4String rangeTypeLabel; 156 G4String rangeTypeLabel;
166 G4String triggerTypeLabel; 157 G4String triggerTypeLabel;
167 G4String decorTriggerTypeLabel; 158 G4String decorTriggerTypeLabel;
168 159
169 if ( data.isARHisto ) 160 if ( data.isARHisto )
170 { 161 {
171 if ( data.isARRec ) 162 if ( data.isARRec )
172 { 163 {
173 fullName += "_arrec"; 164 fullName += "_arrec";
174 rangeTypeLabel = "rec"; 165 rangeTypeLabel = "rec";
175 } 166 }
176 else 167 else
177 { 168 {
178 fullName += "_arreal"; 169 fullName += "_arreal";
179 rangeTypeLabel = "real"; 170 rangeTypeLabel = "real";
180 } 171 }
181 } 172 }
182 173
183 switch ( data.triggerType ) 174 switch ( data.triggerType )
184 { 175 {
185 case CexmcTPT : 176 case CexmcTPT :
186 fullName += "_tpt"; 177 fullName += "_tpt";
187 decorTriggerTypeLabel = " --tpt--"; 178 decorTriggerTypeLabel = " --tpt--";
188 fullTitle += decorTriggerTypeLabel; 179 fullTitle += decorTriggerTypeLabel;
189 triggerTypeLabel = "tpt"; 180 triggerTypeLabel = "tpt";
190 break; 181 break;
191 case CexmcEDT : 182 case CexmcEDT :
192 fullName += "_edt"; 183 fullName += "_edt";
193 decorTriggerTypeLabel = " --edt--"; 184 decorTriggerTypeLabel = " --edt--";
194 fullTitle += decorTriggerTypeLabel; 185 fullTitle += decorTriggerTypeLabel;
195 triggerTypeLabel = "edt"; 186 triggerTypeLabel = "edt";
196 break; 187 break;
197 case CexmcRT : 188 case CexmcRT :
198 fullName += "_rt"; 189 fullName += "_rt";
199 decorTriggerTypeLabel = " --rt--"; 190 decorTriggerTypeLabel = " --rt--";
200 fullTitle += decorTriggerTypeLabel; 191 fullTitle += decorTriggerTypeLabel;
201 triggerTypeLabel = "rt"; 192 triggerTypeLabel = "rt";
202 break; 193 break;
203 default : 194 default :
204 break; 195 break;
205 } 196 }
206 197
207 CexmcHistosMap::iterator found( histos.fi 198 CexmcHistosMap::iterator found( histos.find( data.type ) );
208 << 199
209 if ( found == histos.end() ) 200 if ( found == histos.end() )
210 throw CexmcException( CexmcWeirdExcept 201 throw CexmcException( CexmcWeirdException );
211 202
212 CexmcHistoVector & histoVector( found->se 203 CexmcHistoVector & histoVector( found->second );
213 204
214 if ( data.isARHisto ) 205 if ( data.isARHisto )
215 { 206 {
216 G4bool dirOk( false ); 207 G4bool dirOk( false );
217 208
218 dirOk = gDirectory->Get( fullName ) != << 209 if ( outFile )
>> 210 {
>> 211 dirOk = gDirectory->Get( fullName ) != NULL;
219 212
220 if ( ! dirOk ) << 213 if ( ! dirOk )
221 dirOk = ( gDirectory->mkdir( fullN << 214 dirOk = ( gDirectory->mkdir( fullName, fullTitle ) != NULL );
222 215
223 if ( dirOk ) << 216 if ( dirOk )
224 gDirectory->cd( fullName ); << 217 gDirectory->cd( fullName );
>> 218 }
225 219
226 std::ostringstream histoName; 220 std::ostringstream histoName;
227 std::ostringstream histoTitle; 221 std::ostringstream histoTitle;
228 histoName << data.name << "_r" << aRan 222 histoName << data.name << "_r" << aRange.index + 1 << rangeTypeLabel <<
229 "_" << triggerTypeLabel; 223 "_" << triggerTypeLabel;
230 histoTitle << data.title << " {range " 224 histoTitle << data.title << " {range " << aRange.index + 1 <<
231 rangeTypeLabel << " [" << std: 225 rangeTypeLabel << " [" << std::fixed <<
232 std::setprecision( 4 ) << aRan 226 std::setprecision( 4 ) << aRange.top << ", " <<
233 aRange.bottom << ")}" << decor 227 aRange.bottom << ")}" << decorTriggerTypeLabel;
234 CreateHisto( histoVector, data.impl, h 228 CreateHisto( histoVector, data.impl, histoName.str(), histoTitle.str(),
235 data.axes ); 229 data.axes );
236 230
237 if ( dirOk ) << 231 if ( outFile )
238 gDirectory->cd( ".." ); << 232 {
>> 233 if ( dirOk )
>> 234 gDirectory->cd( ".." );
>> 235 }
239 } 236 }
240 else 237 else
241 { 238 {
242 CreateHisto( histoVector, data.impl, f 239 CreateHisto( histoVector, data.impl, fullName, fullTitle, data.axes );
243 } 240 }
244 } 241 }
245 242
246 243
247 void CexmcHistoManager::CreateHisto( CexmcHis 244 void CexmcHistoManager::CreateHisto( CexmcHistoVector & histoVector,
248 CexmcHistoImpl histoI 245 CexmcHistoImpl histoImpl, const G4String & name,
249 const G4String & titl 246 const G4String & title, const CexmcHistoAxes & axes )
250 { 247 {
251 TH1 * histo( NULL ); 248 TH1 * histo( NULL );
252 249
253 switch ( histoImpl ) 250 switch ( histoImpl )
254 { 251 {
255 case Cexmc_TH1F : 252 case Cexmc_TH1F :
256 histo = new TH1F( name, title, axes.at 253 histo = new TH1F( name, title, axes.at( 0 ).nBins,
257 axes.at( 0 ).nBinsMi 254 axes.at( 0 ).nBinsMin, axes.at( 0 ).nBinsMax );
258 break; 255 break;
259 case Cexmc_TH2F : 256 case Cexmc_TH2F :
260 histo = new TH2F( name, title, axes.at 257 histo = new TH2F( name, title, axes.at( 0 ).nBins,
261 axes.at( 0 ).nBinsMi 258 axes.at( 0 ).nBinsMin, axes.at( 0 ).nBinsMax,
262 axes.at( 1 ).nBins, 259 axes.at( 1 ).nBins, axes.at( 1 ).nBinsMin,
263 axes.at( 1 ).nBinsMa 260 axes.at( 1 ).nBinsMax );
264 break; 261 break;
265 case Cexmc_TH3F : 262 case Cexmc_TH3F :
266 histo = new TH3F( name, title, axes.at 263 histo = new TH3F( name, title, axes.at( 0 ).nBins,
267 axes.at( 0 ).nBinsMi 264 axes.at( 0 ).nBinsMin, axes.at( 0 ).nBinsMax,
268 axes.at( 1 ).nBins, 265 axes.at( 1 ).nBins, axes.at( 1 ).nBinsMin,
269 axes.at( 1 ).nBinsMa 266 axes.at( 1 ).nBinsMax, axes.at( 2 ).nBins,
270 axes.at( 2 ).nBinsMi 267 axes.at( 2 ).nBinsMin, axes.at( 2 ).nBinsMax );
271 break; 268 break;
272 default : 269 default :
273 break; 270 break;
274 } 271 }
275 272
276 if ( histo ) 273 if ( histo )
277 histoVector.push_back( histo ); 274 histoVector.push_back( histo );
278 } 275 }
279 276
280 277
281 void CexmcHistoManager::Initialize( void ) 278 void CexmcHistoManager::Initialize( void )
282 { 279 {
283 if ( isInitialized ) 280 if ( isInitialized )
284 return; 281 return;
285 282
286 CexmcRunManager * runManager( static 283 CexmcRunManager * runManager( static_cast< CexmcRunManager * >(
287 G4 284 G4RunManager::GetRunManager() ) );
288 CexmcPhysicsManager * physicsManager( ru 285 CexmcPhysicsManager * physicsManager( runManager->GetPhysicsManager() );
289 286
290 if ( ! physicsManager ) 287 if ( ! physicsManager )
291 throw CexmcException( CexmcWeirdExcept 288 throw CexmcException( CexmcWeirdException );
292 289
293 CexmcProductionModel * productionModel( p 290 CexmcProductionModel * productionModel( physicsManager->
294 291 GetProductionModel() );
295 292
296 if ( ! productionModel ) 293 if ( ! productionModel )
297 throw CexmcException( CexmcWeirdExcept 294 throw CexmcException( CexmcWeirdException );
298 295
299 G4ParticleDefinition * outputParticle( 296 G4ParticleDefinition * outputParticle(
300 productionMode 297 productionModel->GetOutputParticle() );
301 G4ParticleDefinition * nucleusOutputParti 298 G4ParticleDefinition * nucleusOutputParticle(
302 productionMode 299 productionModel->GetNucleusOutputParticle() );
303 300
304 if ( ! outputParticle || ! nucleusOutputPa 301 if ( ! outputParticle || ! nucleusOutputParticle )
305 throw CexmcException( CexmcIncompleteP 302 throw CexmcException( CexmcIncompleteProductionModel );
306 303
307 opName = outputParticle->GetParticleName() 304 opName = outputParticle->GetParticleName();
308 nopName = nucleusOutputParticle->GetPartic 305 nopName = nucleusOutputParticle->GetParticleName();
309 opMass = outputParticle->GetPDGMass(); 306 opMass = outputParticle->GetPDGMass();
310 nopMass = nucleusOutputParticle->GetPDGMas 307 nopMass = nucleusOutputParticle->GetPDGMass();
311 308
312 G4String title; 309 G4String title;
313 Int_t nBinsX; 310 Int_t nBinsX;
314 Int_t nBinsY; 311 Int_t nBinsY;
315 Double_t nBinsMinX; 312 Double_t nBinsMinX;
316 Double_t nBinsMaxX; 313 Double_t nBinsMaxX;
317 Double_t nBinsMinY; 314 Double_t nBinsMinY;
318 Double_t nBinsMaxY; 315 Double_t nBinsMaxY;
319 CexmcHistoAxes axes; 316 CexmcHistoAxes axes;
320 317
321 if ( runManager->ProjectIsSaved() ) 318 if ( runManager->ProjectIsSaved() )
322 { 319 {
323 G4String projectsDir( runManager->Get 320 G4String projectsDir( runManager->GetProjectsDir() );
324 G4String resultsFile( projectsDir + " 321 G4String resultsFile( projectsDir + "/" + runManager->GetProjectId() +
325 ".root" ); 322 ".root" );
326 outFile = new TFile( resultsFile, "rec 323 outFile = new TFile( resultsFile, "recreate" );
327 } 324 }
328 else <<
329 { <<
330 outFile = new TDirectoryFile( CexmcHis <<
331 CexmcHis <<
332 gDirectory->cd( CexmcHistoDirectoryHan <<
333 } <<
334 <<
335 if ( ! outFile ) <<
336 throw CexmcException( CexmcWeirdExcept <<
337 325
338 const CexmcSetup * setup( static_cast< co 326 const CexmcSetup * setup( static_cast< const CexmcSetup * >(
339 runManager->Ge 327 runManager->GetUserDetectorConstruction() ) );
340 if ( ! setup ) 328 if ( ! setup )
341 throw CexmcException( CexmcWeirdExcept 329 throw CexmcException( CexmcWeirdException );
342 330
343 const G4LogicalVolume * lVolume( setup->G 331 const G4LogicalVolume * lVolume( setup->GetVolume( CexmcSetup::Monitor ) );
344 332
345 if ( ! lVolume ) 333 if ( ! lVolume )
346 throw CexmcException( CexmcIncompatibl 334 throw CexmcException( CexmcIncompatibleGeometry );
347 335
348 nBinsMinX = CexmcHistoBeamMomentumMin; 336 nBinsMinX = CexmcHistoBeamMomentumMin;
349 nBinsMaxX = CexmcHistoBeamMomentumMax; 337 nBinsMaxX = CexmcHistoBeamMomentumMax;
350 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 338 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) /
351 CexmcHistoBeamMomentumReso 339 CexmcHistoBeamMomentumResolution );
352 axes.push_back( CexmcHistoAxisData( nBinsX 340 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
353 AddHisto( CexmcHistoData( CexmcMomentumBP_ 341 AddHisto( CexmcHistoData( CexmcMomentumBP_TPT_Histo, Cexmc_TH1F, false,
354 false, CexmcTPT, "mombp", "Beam mome 342 false, CexmcTPT, "mombp", "Beam momentum at the monitor", axes ) );
355 AddHisto( CexmcHistoData( CexmcMomentumBP_ 343 AddHisto( CexmcHistoData( CexmcMomentumBP_RT_Histo, Cexmc_TH1F, false,
356 false, CexmcRT, "mombp", "Beam momen 344 false, CexmcRT, "mombp", "Beam momentum at the monitor", axes ) );
357 if ( verboseLevel > 0 ) 345 if ( verboseLevel > 0 )
358 { 346 {
359 AddHisto( CexmcHistoData( CexmcMomentu 347 AddHisto( CexmcHistoData( CexmcMomentumIP_TPT_Histo, Cexmc_TH1F, false,
360 false, CexmcTPT, "momip", "Momentu 348 false, CexmcTPT, "momip", "Momentum of the incident particle",
361 axes ) ); 349 axes ) );
362 } 350 }
363 351
364 G4Box * box( dynamic_cast< G4Box * >( lV 352 G4Box * box( dynamic_cast< G4Box * >( lVolume->GetSolid() ) );
365 353
366 if ( ! box ) 354 if ( ! box )
367 throw CexmcException( CexmcIncompatibl 355 throw CexmcException( CexmcIncompatibleGeometry );
368 356
369 G4double width( box->GetXHalfLength() * 2 357 G4double width( box->GetXHalfLength() * 2 );
370 G4double height( box->GetYHalfLength() * 358 G4double height( box->GetYHalfLength() * 2 );
371 G4double halfWidth( width / 2 + CexmcHist 359 G4double halfWidth( width / 2 + CexmcHistoTPSafetyArea );
372 G4double halfHeight( height / 2 + CexmcHi 360 G4double halfHeight( height / 2 + CexmcHistoTPSafetyArea );
373 361
374 nBinsX = Int_t( halfWidth * 2 / CexmcHisto 362 nBinsX = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
375 nBinsY = Int_t( halfHeight * 2 / CexmcHist 363 nBinsY = Int_t( halfHeight * 2 / CexmcHistoTPResolution );
376 axes.clear(); 364 axes.clear();
377 axes.push_back( CexmcHistoAxisData( nBinsX 365 axes.push_back( CexmcHistoAxisData( nBinsX, -halfWidth, halfWidth ) );
378 axes.push_back( CexmcHistoAxisData( nBinsY 366 axes.push_back( CexmcHistoAxisData( nBinsY, -halfHeight, halfHeight ) );
379 AddHisto( CexmcHistoData( CexmcTPInMonitor 367 AddHisto( CexmcHistoData( CexmcTPInMonitor_TPT_Histo, Cexmc_TH2F, false,
380 false, CexmcTPT, "tpmon", "Track point 368 false, CexmcTPT, "tpmon", "Track points (mon)", axes ) );
381 369
382 lVolume = setup->GetVolume( CexmcSetup::Ta 370 lVolume = setup->GetVolume( CexmcSetup::Target );
383 G4Tubs * tube( dynamic_cast< G4Tubs * >( 371 G4Tubs * tube( dynamic_cast< G4Tubs * >( lVolume->GetSolid() ) );
384 372
385 if ( ! tube ) 373 if ( ! tube )
386 throw CexmcException( CexmcIncompatibl 374 throw CexmcException( CexmcIncompatibleGeometry );
387 375
388 G4double radius( tube->GetOuterRadius() ) 376 G4double radius( tube->GetOuterRadius() );
389 height = tube->GetZHalfLength() * 2; 377 height = tube->GetZHalfLength() * 2;
390 halfWidth = radius + CexmcHistoTPSafetyAre 378 halfWidth = radius + CexmcHistoTPSafetyArea;
391 halfHeight = height / 2 + CexmcHistoTPSafe 379 halfHeight = height / 2 + CexmcHistoTPSafetyArea;
392 380
393 nBinsX = Int_t( halfWidth * 2 / CexmcHisto 381 nBinsX = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
394 nBinsY = Int_t( halfWidth * 2 / CexmcHisto 382 nBinsY = Int_t( halfWidth * 2 / CexmcHistoTPResolution );
395 Int_t nBinsZ( Int_t( halfHeight * 2 / Cex 383 Int_t nBinsZ( Int_t( halfHeight * 2 / CexmcHistoTPResolution ) );
396 axes.clear(); 384 axes.clear();
397 axes.push_back( CexmcHistoAxisData( nBinsX 385 axes.push_back( CexmcHistoAxisData( nBinsX, -halfWidth, halfWidth ) );
398 axes.push_back( CexmcHistoAxisData( nBinsY 386 axes.push_back( CexmcHistoAxisData( nBinsY, -halfWidth, halfWidth ) );
399 axes.push_back( CexmcHistoAxisData( nBinsZ 387 axes.push_back( CexmcHistoAxisData( nBinsZ, -halfHeight, halfHeight ) );
400 AddHisto( CexmcHistoData( CexmcTPInTarget_ 388 AddHisto( CexmcHistoData( CexmcTPInTarget_TPT_Histo, Cexmc_TH3F, false,
401 false, CexmcTPT, "tptar", "Track point 389 false, CexmcTPT, "tptar", "Track points (tar)", axes ) );
402 AddHisto( CexmcHistoData( CexmcTPInTarget_ 390 AddHisto( CexmcHistoData( CexmcTPInTarget_RT_Histo, Cexmc_TH3F, false,
403 false, CexmcRT, "tptar", "Track points 391 false, CexmcRT, "tptar", "Track points (tar)", axes ) );
404 392
405 title = "Reconstructed masses (" + nopName 393 title = "Reconstructed masses (" + nopName + " vs. " + opName + ")";
406 nBinsMinX = opMass / 2; 394 nBinsMinX = opMass / 2;
407 nBinsMaxX = opMass + opMass / 2; 395 nBinsMaxX = opMass + opMass / 2;
408 nBinsMinY = nopMass / 2; 396 nBinsMinY = nopMass / 2;
409 nBinsMaxY = nopMass + nopMass / 2; 397 nBinsMaxY = nopMass + nopMass / 2;
410 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 398 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoMassResolution );
411 nBinsY = Int_t( ( nBinsMaxY - nBinsMinY ) 399 nBinsY = Int_t( ( nBinsMaxY - nBinsMinY ) / CexmcHistoMassResolution );
412 axes.clear(); 400 axes.clear();
413 axes.push_back( CexmcHistoAxisData( nBinsX 401 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
414 axes.push_back( CexmcHistoAxisData( nBinsY 402 axes.push_back( CexmcHistoAxisData( nBinsY, nBinsMinY, nBinsMaxY ) );
415 AddHisto( CexmcHistoData( CexmcRecMasses_E 403 AddHisto( CexmcHistoData( CexmcRecMasses_EDT_Histo, Cexmc_TH2F, false,
416 false, CexmcEDT, "recmasses", title, a 404 false, CexmcEDT, "recmasses", title, axes ) );
417 AddHisto( CexmcHistoData( CexmcRecMasses_R 405 AddHisto( CexmcHistoData( CexmcRecMasses_RT_Histo, Cexmc_TH2F, false,
418 false, CexmcRT, "recmasses", title, ax 406 false, CexmcRT, "recmasses", title, axes ) );
419 407
420 nBinsMinX = 0.; 408 nBinsMinX = 0.;
421 nBinsMaxX = CexmcHistoEnergyMax; 409 nBinsMaxX = CexmcHistoEnergyMax;
422 nBinsMinY = 0.; 410 nBinsMinY = 0.;
423 nBinsMaxY = CexmcHistoEnergyMax; 411 nBinsMaxY = CexmcHistoEnergyMax;
424 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) 412 nBinsX = Int_t( ( nBinsMaxX - nBinsMinX ) / CexmcHistoEnergyResolution );
425 nBinsY = Int_t( ( nBinsMaxY - nBinsMinY ) 413 nBinsY = Int_t( ( nBinsMaxY - nBinsMinY ) / CexmcHistoEnergyResolution );
426 axes.clear(); 414 axes.clear();
427 axes.push_back( CexmcHistoAxisData( nBinsX 415 axes.push_back( CexmcHistoAxisData( nBinsX, nBinsMinX, nBinsMaxX ) );
428 axes.push_back( CexmcHistoAxisData( nBinsY 416 axes.push_back( CexmcHistoAxisData( nBinsY, nBinsMinY, nBinsMaxY ) );
429 AddHisto( CexmcHistoData( CexmcAbsorbedEne 417 AddHisto( CexmcHistoData( CexmcAbsorbedEnergy_EDT_Histo, Cexmc_TH2F, false,
430 false, CexmcEDT, "ae", "Absorbed energ 418 false, CexmcEDT, "ae", "Absorbed energy (rc vs. lc)", axes ) );
431 AddHisto( CexmcHistoData( CexmcAbsorbedEne 419 AddHisto( CexmcHistoData( CexmcAbsorbedEnergy_RT_Histo, Cexmc_TH2F, false,
432 false, CexmcRT, "ae", "Absorbed energy 420 false, CexmcRT, "ae", "Absorbed energy (rc vs. lc)", axes ) );
433 421
434 SetupARHistos( runManager->GetPhysicsManag 422 SetupARHistos( runManager->GetPhysicsManager()->GetProductionModel()->
435 GetAngularRanges() ); 423 GetAngularRanges() );
436 424
437 isInitialized = true; 425 isInitialized = true;
438 } 426 }
439 427
440 428
441 void CexmcHistoManager::SetupARHistos( const 429 void CexmcHistoManager::SetupARHistos( const CexmcAngularRangeList & aRanges )
442 { 430 {
443 TIter objs( gDirectory->GetList() ); 431 TIter objs( gDirectory->GetList() );
444 TObject * obj( NULL ); 432 TObject * obj( NULL );
445 433
446 while ( ( obj = ( TObject * )objs() ) ) 434 while ( ( obj = ( TObject * )objs() ) )
447 { 435 {
448 TString name( obj->GetName() ); 436 TString name( obj->GetName() );
449 437
450 if ( obj->IsFolder() && ( name.Contain << 438 if ( outFile )
451 name.Contain <<
452 { 439 {
453 gDirectory->cd( name ); << 440 if ( obj->IsFolder() && ( name.Contains( TString( "_arreal_" ) ) ||
454 gDirectory->DeleteAll(); << 441 name.Contains( TString( "_arrec_" ) ) ) )
455 gDirectory->cd( ".." ); << 442 {
>> 443 gDirectory->cd( name );
>> 444 gDirectory->DeleteAll();
>> 445 gDirectory->cd( ".." );
>> 446 }
>> 447 }
>> 448 else
>> 449 {
>> 450 if ( name.Contains( TRegexp( "_r[0-9]+" ) ) )
>> 451 {
>> 452 gDirectory->Remove( obj );
>> 453 }
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 <<
739 gDirectory->ls(); 736 gDirectory->ls();
740 } 737 }
741 738
742 739
743 void CexmcHistoManager::Print( const G4String 740 void CexmcHistoManager::Print( const G4String & value )
744 { 741 {
745 TObject * histo( gDirectory->FindObjectAn 742 TObject * histo( gDirectory->FindObjectAny( value.c_str() ) );
746 743
747 if ( ! histo ) 744 if ( ! histo )
748 { 745 {
749 G4cout << "Histogram '" << value << "' 746 G4cout << "Histogram '" << value << "' was not found" << G4endl;
750 return; 747 return;
751 } 748 }
752 749
753 /* BEWARE: histo will be printed on stdout <<
754 histo->Print( "range" ); 750 histo->Print( "range" );
755 } 751 }
756 752
757 753
758 #ifdef CEXMC_USE_ROOTQT 754 #ifdef CEXMC_USE_ROOTQT
759 755
760 void CexmcHistoManager::Draw( const G4String 756 void CexmcHistoManager::Draw( const G4String & histoName,
761 const G4String 757 const G4String & histoDrawOptions )
762 { 758 {
763 if ( ! areLiveHistogramsEnabled ) << 759 CexmcRunManager * runManager( static_cast< CexmcRunManager * >(
>> 760 G4RunManager::GetRunManager() ) );
>> 761 if ( ! runManager->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 } <<
791 <<
792 <<
793 void CexmcHistoManager::EnableLiveHistograms( <<
794 <<
795 { <<
796 areLiveHistogramsEnabled = on; <<
797 <<
798 if ( ! on || ! isInitialized ) <<
799 return; <<
800 <<
801 G4UIQt * qtSession( dynamic_cast< G4UIQt <<
802 <<
803 if ( ! qtSession ) <<
804 return; <<
805 <<
806 if ( ! histoMenuHandle.empty() && ! isHist <<
807 { <<
808 qtSession->AddMenu( histoMenuHandle, h <<
809 BuildMenuTree( qtSession, histoMenuHan <<
810 isHistoMenuInitialized = true; <<
811 } <<
812 } <<
813 <<
814 <<
815 void CexmcHistoManager::BuildMenuTree( G4UIQt <<
816 const <<
817 { <<
818 TIter objs( ls ); <<
819 TObject * obj( NULL ); <<
820 <<
821 while ( ( obj = ( TObject * )objs() ) ) <<
822 { <<
823 G4String name( obj->GetName() ); <<
824 G4String title( obj->GetTitle() ); <<
825 <<
826 if ( obj->IsFolder() ) <<
827 { <<
828 AddSubmenu( session, menu, name, t <<
829 BuildMenuTree( session, name, ( ( <<
830 } <<
831 else <<
832 { <<
833 G4String options( name ); <<
834 <<
835 do <<
836 { <<
837 if ( obj->InheritsFrom( TH3::C <<
838 ! drawOptions3D.empty() ) <<
839 { <<
840 title = G4String( "3: " ) <<
841 options += G4String( " " ) <<
842 break; <<
843 } <<
844 if ( obj->InheritsFrom( TH2::C <<
845 ! drawOptions2D.empty() ) <<
846 { <<
847 title = G4String( "2: " ) <<
848 options += G4String( " " ) <<
849 break; <<
850 } <<
851 if ( obj->InheritsFrom( TH1::C <<
852 ! drawOptions1D.empty() ) <<
853 { <<
854 options += G4String( " " ) <<
855 break; <<
856 } <<
857 } while ( false ); <<
858 <<
859 G4String cmd( CexmcMessenger::his <<
860 session->AddButton( menu, title.c_ <<
861 } <<
862 } <<
863 } <<
864 <<
865 <<
866 void CexmcHistoManager::AddSubmenu( G4UIQt * <<
867 const G4S <<
868 const G4S <<
869 const G4S <<
870 { <<
871 QMenu * menu( new QMenu( label.c_str() ) ); <<
872 QMenu * parentMenu( ( QMenu * )session->Get <<
873 <<
874 parentMenu->addMenu( menu ); <<
875 session->AddInteractor( name, ( G4Interactor <<
876 } 788 }
877 789
878 #endif 790 #endif
879 791
880 #endif 792 #endif
881 793
882 794