source: trunk/source/run/src/G4VUserPhysicsList.cc

Last change on this file was 1337, checked in by garnier, 14 years ago

tag geant4.9.4 beta 1 + modifs locales
File size: 29.0 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer                                           *
4// *                                                                  *
5// * The  Geant4 software  is  copyright of the Copyright Holders  of *
6// * the Geant4 Collaboration.  It is provided  under  the terms  and *
7// * conditions of the Geant4 Software License,  included in the file *
8// * LICENSE and available at  http://cern.ch/geant4/license .  These *
9// * include a list of copyright holders.                             *
10// *                                                                  *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work  make  any representation or  warranty, express or implied, *
14// * regarding  this  software system or assume any liability for its *
15// * use.  Please see the license in the file  LICENSE  and URL above *
16// * for the full disclaimer and the limitation of liability.         *
17// *                                                                  *
18// * This  code  implementation is the result of  the  scientific and *
19// * technical work of the GEANT4 collaboration.                      *
20// * By using,  copying,  modifying or  distributing the software (or *
21// * any work based  on the software)  you  agree  to acknowledge its *
22// * use  in  resulting  scientific  publications,  and indicate your *
23// * acceptance of all terms of the Geant4 Software license.          *
24// ********************************************************************
25//
26//
27// $Id: G4VUserPhysicsList.cc,v 1.71 2009/08/09 14:31:46 kurasige Exp $
28// GEANT4 tag $Name: geant4-09-04-beta-01 $
29//
30//
31// ------------------------------------------------------------
32//      GEANT 4 class header file
33//
34// ------------------------------------------------------------
35//      History
36//       first version                    09 Jan 1998 by H.Kurashige
37//       Added SetEnergyRange             18 Jun 1998 by H.Kurashige
38//       Change for short lived particles 27 Jun 1998 by H.Kurashige
39//       G4BestUnit on output             12 nov 1998 by M.Maire
40//       Added RemoveProcessManager        9 Feb 1999 by H.Kurashige
41//       Fixed RemoveProcessManager       15 Apr 1999 by H.Kurashige
42//       Removed ConstructAllParticles()  15 Apr 1999 by H.Kurashige
43//       Modified for CUTS per REGION     10 Oct 2002 by H.Kurashige
44//       Check if particle IsShortLived   18 Jun 2003 by V.Ivanchenko
45//       Modify PreparePhysicsList        18 Jan 2006 by H.Kurashige
46// ------------------------------------------------------------
47
48#include "globals.hh"
49#include "G4VUserPhysicsList.hh"
50#include "G4ParticleDefinition.hh"
51#include "G4ProcessManager.hh"
52#include "G4ParticleTable.hh"
53#include "G4ProductionCutsTable.hh"
54#include "G4Material.hh"
55#include "G4UserPhysicsListMessenger.hh"
56#include "G4UImanager.hh"
57#include "G4UnitsTable.hh"
58#include "G4RegionStore.hh"
59#include "G4Region.hh"
60#include "G4ProductionCutsTable.hh"
61#include "G4ProductionCuts.hh"
62#include "G4MaterialCutsCouple.hh"
63
64#include "G4ios.hh"
65#include <iomanip>
66#include <fstream>
67
68////////////////////////////////////////////////////////
69G4VUserPhysicsList::G4VUserPhysicsList()
70                 :  fDisableCheckParticleList(false),
71                    verboseLevel(1),
72                    fRetrievePhysicsTable(false),
73                    fStoredInAscii(true),
74                    fIsCheckedForRetrievePhysicsTable(false),
75                    fIsRestoredCutValues(false),
76                    directoryPhysicsTable("."),
77                    fDisplayThreshold(0),
78                    fIsPhysicsTableBuilt(false),
79                    useCoupledTransportation(false)
80{
81  // default cut value  (1.0mm)
82  defaultCutValue = 1.0*mm;
83
84  // pointer to the particle table
85  theParticleTable = G4ParticleTable::GetParticleTable();
86  theParticleIterator = theParticleTable->GetIterator();
87
88  // pointer to the cuts table
89  fCutsTable =  G4ProductionCutsTable::GetProductionCutsTable();
90
91  // set energy range for SetCut calcuration
92  fCutsTable->SetEnergyRange(0.99*keV, 100*TeV);
93
94  // UI Messenger
95  theMessenger = new G4UserPhysicsListMessenger(this);
96}
97
98////////////////////////////////////////////////////////
99G4VUserPhysicsList::~G4VUserPhysicsList()
100{
101  if (theMessenger != 0) {
102    delete theMessenger;
103    theMessenger = 0;
104  }
105  RemoveProcessManager();
106
107  // invoke DeleteAllParticle
108  theParticleTable->DeleteAllParticles();
109
110}
111
112////////////////////////////////////////////////////////
113void G4VUserPhysicsList::SetVerboseLevel(G4int value)
114{
115  verboseLevel = value;
116  // set verboseLevel for G4ProductionCutsTable same as one for G4VUserPhysicsList:
117  fCutsTable->SetVerboseLevel(value);
118
119#ifdef G4VERBOSE
120  if (verboseLevel >1){
121    G4cout << "G4VUserPhysicsList::SetVerboseLevel  :";
122    G4cout << " Verbose level is set to " << verboseLevel << G4endl;
123  }
124#endif
125}
126
127////////////////////////////////////////////////////////
128void G4VUserPhysicsList::AddProcessManager(G4ParticleDefinition* newParticle,
129                                           G4ProcessManager*     newManager)
130{
131  if (newParticle == 0) return;
132  if (newParticle->GetProcessManager() != 0) {
133#ifdef G4VERBOSE
134    if (verboseLevel >1){
135      G4cout << "G4VUserPhysicsList::AddProcessManager: ";
136      G4cout  << newParticle->GetParticleName();
137      G4cout << " already has ProcessManager " << G4endl;
138    }
139#endif
140    return;
141  }
142
143  // create new process manager if newManager  == 0
144  if (newManager  == 0){
145    // Add ProcessManager
146    if (newParticle->GetParticleType() == "nucleus") {
147      // Create a copy of the process manager of "GenericIon" in case of "nucleus"
148      G4ParticleDefinition* genericIon =
149           (G4ParticleTable::GetParticleTable())->FindParticle("GenericIon");
150
151      if (genericIon != 0) {
152        G4ProcessManager* ionMan = genericIon->GetProcessManager();
153        if (ionMan != 0) {
154          newManager = new G4ProcessManager(*ionMan);
155        } else {
156          // no process manager has been registered yet
157          newManager = new G4ProcessManager(newParticle);
158          G4Exception("G4VUserPhysicsList::AddProcessManager","Error in GenericIon",
159                RunMustBeAborted,"GenericIon has no ProcessMamanger"); 
160        }
161      } else {
162        // "GenericIon" does not exist
163        newManager = new G4ProcessManager(newParticle);
164        G4Exception("G4VUserPhysicsList::AddProcessManager","No GenericIon",
165                    RunMustBeAborted,"GenericIon does not exist");     
166      }
167
168    } else {
169      // create process manager for particles other than "nucleus"
170      newManager = new G4ProcessManager(newParticle);
171    }
172  }
173
174  // set particle type
175  newManager->SetParticleType(newParticle);
176
177  // add the process manager
178  newParticle->SetProcessManager(newManager);
179
180#ifdef G4VERBOSE
181 if (verboseLevel >2){
182    G4cout << "G4VUserPhysicsList::AddProcessManager: ";
183    G4cout  << "adds ProcessManager to ";
184    G4cout  << newParticle->GetParticleName() << G4endl;
185    newManager->DumpInfo();
186  }
187#endif
188  if ( fIsPhysicsTableBuilt
189       && (newParticle->GetParticleType() == "nucleus")) {
190    PreparePhysicsTable(newParticle);
191    BuildPhysicsTable(newParticle);
192  }
193}
194
195
196////////////////////////////////////////////////////////
197void G4VUserPhysicsList::CheckParticleList()
198{
199
200  // skip if fDisableCheckParticleList is set 
201  if (fDisableCheckParticleList) return;
202
203  bool isElectron = false;
204  bool isPositron = false;
205  bool isGamma    = false;
206  bool isProton   = false;
207  bool isGenericIon = false;
208  bool isAnyIon   = false;
209  bool isAnyChargedBaryon   = false;
210  bool isEmProc   = false;
211
212  // loop over all particles in G4ParticleTable
213  theParticleIterator->reset();
214  while( (*theParticleIterator)() ){
215    G4ParticleDefinition* particle = theParticleIterator->value();
216    G4String name = particle->GetParticleName();
217    // check if any EM process exists
218    if (!isEmProc) {
219      G4ProcessVector* list = particle->GetProcessManager()->GetProcessList();
220      for (int idx=0; idx<list->size(); idx++){
221        isEmProc = ((*list)[idx])->GetProcessType() == fElectromagnetic;
222        if (isEmProc) break;
223      }
224    }
225   
226    if      ( name == "e-") isElectron = true; 
227    else if ( name == "e+") isPositron = true; 
228    else if ( name == "gamma") isGamma = true; 
229    else if ( name == "GenericIon") isGenericIon = true; 
230    else if ( name == "proton") isProton = true; 
231    else if ( particle->GetParticleType() == "nucleus") isAnyIon = true;
232    else if ( particle->GetParticleType() == "baryon") {
233       if ( particle->GetPDGCharge() != 0.0 ) isAnyChargedBaryon = true;
234    }
235  }
236
237  if (!isEmProc) return;
238
239  // RULE 1
240  //  e+, e- and gamma should exist
241  //   if one of them exist
242  bool isEmBasic =  isElectron || isPositron || isGamma;
243  bool isMissingEmBasic =  !isElectron || !isPositron || !isGamma;
244  if (isEmBasic && isMissingEmBasic) {
245    G4String missingName="";
246    if (!isElectron) missingName += "e- ";
247    if (!isPositron) missingName += "e+ ";
248    if (!isGamma) missingName += "gamma ";
249
250#ifdef G4VERBOSE
251    if (verboseLevel >0){
252      G4cout << "G4VUserPhysicsList::CheckParticleList: ";
253      G4cout << missingName << " do not exist " << G4endl; 
254      G4cout << " These particle are necessary for basic EM processes" << G4endl;
255    }
256#endif
257    missingName += " should be created ";
258    G4Exception("G4VUserPhysicsList::CheckParticleList","Missing EM basic particle",
259                FatalException, missingName);   
260  }
261
262  // RULE 2
263  //  proton should exist
264  //   if any other charged baryon  exist
265  if (!isProton && isAnyChargedBaryon) {
266    G4String missingName="proton ";
267
268#ifdef G4VERBOSE
269    if (verboseLevel >0){
270      G4cout << "G4VUserPhysicsList::CheckParticleList: ";
271      G4cout << missingName << " does not exist "<< G4endl; 
272      G4cout << " Proton is necessary for EM baryon processes" << G4endl;
273    }
274#endif
275    missingName += " should be created ";
276    G4Exception("G4VUserPhysicsList::CheckParticleList","Missing Proton",
277                FatalException, missingName);   
278  }
279   
280  // RULE 3
281  //  GenericIonn should exist
282  //   if any other ion  exist
283  if (!isGenericIon && isAnyIon) {
284    G4String missingName="GenericIon ";
285
286#ifdef G4VERBOSE
287    if (verboseLevel >0){
288      G4cout << "G4VUserPhysicsList::CheckParticleList: ";
289      G4cout << missingName << " does not exist "<< G4endl; 
290      G4cout << " GenericIon should be created if any ion is necessary" << G4endl;
291    }
292#endif
293    missingName += " should be created ";
294    G4Exception("G4VUserPhysicsList::CheckParticleList","Missing GenericIon",
295                FatalException, missingName);   
296  }
297     
298}
299
300////////////////////////////////////////////////////////
301void G4VUserPhysicsList::InitializeProcessManager()
302{
303  // loop over all particles in G4ParticleTable
304  theParticleIterator->reset();
305  while( (*theParticleIterator)() ){
306    G4ParticleDefinition* particle = theParticleIterator->value();
307    G4ProcessManager* pmanager = particle->GetProcessManager();
308    if  (pmanager==0) {
309      // create process manager if the particle has no its one
310      pmanager = new G4ProcessManager(particle);
311      particle->SetProcessManager(pmanager);
312    }
313  }
314}
315
316/////////////////////////////////////////////////////////
317void G4VUserPhysicsList::RemoveProcessManager()
318{
319  // loop over all particles in G4ParticleTable
320  theParticleIterator->reset();
321  while( (*theParticleIterator)() ){
322    G4ParticleDefinition* particle = theParticleIterator->value();
323    G4ProcessManager* pmanager = particle->GetProcessManager();
324    if  (pmanager!=0) delete pmanager;
325    particle->SetProcessManager(0);
326#ifdef G4VERBOSE
327    if (verboseLevel >2){
328      G4cout << "G4VUserPhysicsList::RemoveProcessManager: ";
329      G4cout  << "remove ProcessManager from ";
330      G4cout  << particle->GetParticleName() << G4endl;
331    }
332#endif
333  }
334}
335
336
337////////////////////////////////////////////////////////
338#include "G4Transportation.hh"
339#include "G4CoupledTransportation.hh"
340#include "G4RunManagerKernel.hh"
341#include "G4ScoringManager.hh"
342
343void G4VUserPhysicsList::AddTransportation()
344{
345  G4int verboseLevelTransport = 0;
346  G4VProcess* theTransportationProcess = 0;
347  G4int nParaWorld = G4RunManagerKernel::GetRunManagerKernel()->GetNumberOfParallelWorld();
348
349  if(nParaWorld || useCoupledTransportation || G4ScoringManager::GetScoringManagerIfExist())
350    {
351      theTransportationProcess = new G4CoupledTransportation(verboseLevelTransport);
352      G4cout << "#############################################################################" << G4endl
353             << " G4VUserPhysicsList::AddTransportation() --- G4CoupledTransportation is used " << G4endl
354             << "#############################################################################" << G4endl;
355    }
356  else
357    { theTransportationProcess = new G4Transportation(verboseLevelTransport); }
358 
359#ifdef G4VERBOSE
360    if (verboseLevel >2){
361      G4cout << "G4VUserPhysicsList::AddTransportation()  "<< G4endl;
362    }
363#endif
364
365  // loop over all particles in G4ParticleTable
366  theParticleIterator->reset();
367  while( (*theParticleIterator)() ){
368    G4ParticleDefinition* particle = theParticleIterator->value();
369    G4ProcessManager* pmanager = particle->GetProcessManager();
370    if (!particle->IsShortLived()) {
371      // Add transportation process for all particles other than  "shortlived"
372      if ( pmanager == 0) {
373        // Error !! no process manager
374        G4String particleName = particle->GetParticleName();
375        G4Exception("G4VUserPhysicsList::AddTransportation","No process manager",
376                    FatalException, particleName );
377      } else {
378        // add transportation with ordering = ( -1, "first", "first" )
379        pmanager ->AddProcess(theTransportationProcess);
380        pmanager ->SetProcessOrderingToFirst(theTransportationProcess, idxAlongStep);
381        pmanager ->SetProcessOrderingToFirst(theTransportationProcess, idxPostStep);
382      }
383    } else {
384      // shortlived particle case
385     // Add transportation process for all particles other than  "shortlived"
386      if ( pmanager == 0) {
387        // Error !! no process manager
388        G4String particleName = particle->GetParticleName();
389        G4Exception("G4VUserPhysicsList::AddTransportation","No process manager",
390                    FatalException, particleName );
391      } else {
392        // add transportation with ordering = ( -1, "first", "first" )
393        pmanager ->AddProcess(theTransportationProcess);
394        pmanager ->SetProcessOrderingToFirst(theTransportationProcess, idxAlongStep);
395        pmanager ->SetProcessOrderingToFirst(theTransportationProcess, idxPostStep);
396      }
397
398    }
399  }
400}
401
402
403////////////////////////////////////////////////////////
404void G4VUserPhysicsList::SetDefaultCutValue(G4double value)
405{
406   if (value<=0.0) {
407#ifdef G4VERBOSE
408     if (verboseLevel >0){
409       G4cout << "G4VUserPhysicsList::SetDefaultCutValue: negative cut values";
410       G4cout << "  :" << value/mm << "[mm]" << G4endl;
411     }
412#endif
413   } else {
414#ifdef G4VERBOSE
415     if (verboseLevel >1){
416       G4cout << "G4VUserPhysicsList::SetDefaultCutValue:";
417       G4cout << "default cut value is changed to   :" ;
418       G4cout << value/mm << "[mm]" << G4endl;
419     }
420#endif
421     defaultCutValue = value;
422
423   }
424}
425
426
427////////////////////////////////////////////////////////
428void G4VUserPhysicsList::SetCutValue(G4double aCut, const G4String& name)
429{
430  G4ParticleDefinition* particle = theParticleTable->FindParticle(name);
431  if (particle != 0){
432    if (!particle->IsShortLived()) {
433      //set cut value
434      SetParticleCuts( aCut ,particle );
435    }
436  }
437}
438
439////////////////////////////////////////////////////////
440void G4VUserPhysicsList::SetCutValue
441(G4double aCut, const G4String& pname, const G4String& rname)
442{
443  G4ParticleDefinition* particle = theParticleTable->FindParticle(pname);
444  G4Region* region = G4RegionStore::GetInstance()->GetRegion(rname);
445  if (particle != 0 && region != 0){
446    if (!particle->IsShortLived()) {
447      //set cut value
448      SetParticleCuts( aCut ,particle, region );
449    }
450  }
451}
452
453
454
455////////////////////////////////////////////////////////
456void G4VUserPhysicsList::SetCutsWithDefault()
457{
458  // default cut value
459  G4double cut = defaultCutValue;
460
461#ifdef G4VERBOSE
462  if (verboseLevel >1){
463    G4cout << "G4VUserPhysicsList::SetCutsWithDefault:";
464    G4cout << "CutLength : " << cut/mm << " (mm)" << G4endl;
465  }
466#endif
467
468  // set cut values for gamma at first and for e- and e+
469  SetCutValue(cut, "gamma");
470  SetCutValue(cut, "e-");
471  SetCutValue(cut, "e+");
472  SetCutValue(cut, "proton");
473
474  // dump Cut values if verboseLevel==3
475  if (verboseLevel>2) {
476    DumpCutValuesTable();
477  }
478}
479
480
481////////////////////////////////////////////////////////
482void G4VUserPhysicsList::SetCutsForRegion(G4double aCut, const G4String& rname)
483{
484  // set cut values for gamma at first and for e- and e+
485  SetCutValue(aCut, "gamma", rname);
486  SetCutValue(aCut, "e-", rname);
487  SetCutValue(aCut, "e+", rname);
488  SetCutValue(aCut, "proton", rname);
489}
490
491
492
493////////////////////////////////////////////////////////
494////////////////////////////////////////////////////////
495void G4VUserPhysicsList::SetParticleCuts( G4double cut, G4ParticleDefinition* particle, G4Region* region)
496{
497  if(!region) region = (*(G4RegionStore::GetInstance()))[0];
498  G4ProductionCuts* pcuts = region->GetProductionCuts();
499  pcuts->SetProductionCut(cut,particle);
500}
501
502///////////////////////////////////////////////////////////////
503void G4VUserPhysicsList::BuildPhysicsTable()
504{
505  //Prepare Physics table for all particles
506  theParticleIterator->reset();
507  while( (*theParticleIterator)() ){
508    G4ParticleDefinition* particle = theParticleIterator->value();
509    PreparePhysicsTable(particle); 
510  }
511
512  // ask processes to prepare physics table
513  if (fRetrievePhysicsTable) {
514    fIsRestoredCutValues = fCutsTable->RetrieveCutsTable(directoryPhysicsTable, fStoredInAscii);
515    // check if retrieve Cut Table successfully
516    if (!fIsRestoredCutValues) {
517#ifdef G4VERBOSE
518      if (verboseLevel>0){
519        G4cout << "G4VUserPhysicsList::BuildPhysicsTable";
520        G4cout << " Retrieve Cut Table failed !!" << G4endl;
521      }
522#endif 
523      G4Exception("G4VUserPhysicsList::BuildPhysicsTable","Fail to Retreive",
524                  RunMustBeAborted,"Production Cut Table can not be retreived");
525    } else {
526#ifdef G4VERBOSE
527      if (verboseLevel>2){
528        G4cout << "G4VUserPhysicsList::BuildPhysicsTable";
529        G4cout << "  Retrieve Cut Table successfully " << G4endl;
530      }
531#endif
532    }     
533  } else {
534#ifdef G4VERBOSE
535    if (verboseLevel>2){
536      G4cout << "G4VUserPhysicsList::BuildPhysicsTable";
537      G4cout << " does not retrieve Cut Table but calculate " << G4endl;
538    } 
539#endif     
540  }
541
542  // Sets a value to particle
543  // set cut values for gamma at first and for e- and e+
544  G4String particleName;
545  G4ParticleDefinition* GammaP = theParticleTable->FindParticle("gamma");
546  if(GammaP) BuildPhysicsTable(GammaP);
547  G4ParticleDefinition* EMinusP = theParticleTable->FindParticle("e-");
548  if(EMinusP) BuildPhysicsTable(EMinusP);
549  G4ParticleDefinition* EPlusP = theParticleTable->FindParticle("e+");
550  if(EPlusP) BuildPhysicsTable(EPlusP);
551  G4ParticleDefinition* ProtonP = theParticleTable->FindParticle("proton");
552  if(ProtonP) BuildPhysicsTable(ProtonP);
553
554
555  theParticleIterator->reset();
556  while( (*theParticleIterator)() ){
557    G4ParticleDefinition* particle = theParticleIterator->value();
558    if( particle!=GammaP && 
559        particle!=EMinusP && 
560        particle!=EPlusP && 
561        particle!=ProtonP   ){
562      BuildPhysicsTable(particle); 
563    }
564  }
565
566  // Set flag
567  fIsPhysicsTableBuilt = true;
568
569}
570///////////////////////////////////////////////////////////////
571void G4VUserPhysicsList::BuildPhysicsTable(G4ParticleDefinition* particle)
572{
573  if (fRetrievePhysicsTable) {
574    if ( !fIsRestoredCutValues){
575      // fail to retreive cut tables
576#ifdef G4VERBOSE
577      if (verboseLevel>0){
578        G4cout << "G4VUserPhysicsList::BuildPhysicsTable  ";
579        G4cout << "Physics table can not be retreived and will be calculated "<< G4endl;
580      }
581#endif
582      fRetrievePhysicsTable = false; 
583
584    } else {
585#ifdef G4VERBOSE
586      if (verboseLevel>2){
587        G4cout << "G4VUserPhysicsList::BuildPhysicsTable  ";
588        G4cout << " Retrieve Physics Table for ";
589        G4cout << particle->GetParticleName() << G4endl;
590      }
591#endif
592      //  Retrieve PhysicsTable from files for proccesses
593      RetrievePhysicsTable(particle, directoryPhysicsTable, fStoredInAscii);
594    }
595  }
596
597#ifdef G4VERBOSE
598  if (verboseLevel>2){
599    G4cout << "G4VUserPhysicsList::BuildPhysicsTable  ";
600    G4cout << "Calculate Physics Table for " << particle->GetParticleName() << G4endl;
601  }
602#endif
603  // Rebuild the physics tables for every process for this particle type
604  // if particle is not ShortLived
605  if(!particle->IsShortLived()) {
606    G4ProcessManager* pManager =  particle->GetProcessManager();
607    if (!pManager) {
608      G4cerr << "G4VUserPhysicsList::BuildPhysicsTable  : No Process Manager for " 
609             << particle->GetParticleName() <<G4endl;
610      G4cerr << particle->GetParticleName() << " should be created in your PhysicsList" <<G4endl;
611      G4Exception("G4VUserPhysicsList::BuildPhysicsTable","No process manager",
612                    FatalException,  particle->GetParticleName() );
613    }
614    G4ProcessVector* pVector = pManager->GetProcessList();
615    if (!pVector) {
616      G4cerr << "G4VUserPhysicsList::BuildPhysicsTable  : No Process Vector for " 
617             << particle->GetParticleName() <<G4endl;
618      G4cerr << particle->GetParticleName() << " should be created in your PhysicsList" <<G4endl;
619      G4Exception("G4VUserPhysicsList::BuildPhysicsTable","No process Vector",
620                    FatalException,  particle->GetParticleName() );
621    }
622    for (G4int j=0; j < pVector->size(); ++j) {
623      (*pVector)[j]->BuildPhysicsTable(*particle);
624    }
625  }
626}
627
628///////////////////////////////////////////////////////////////
629void G4VUserPhysicsList::PreparePhysicsTable(G4ParticleDefinition* particle)
630{
631  // Prepare the physics tables for every process for this particle type
632  // if particle is not ShortLived
633  if(!particle->IsShortLived()) {
634    G4ProcessManager* pManager =  particle->GetProcessManager();
635    if (!pManager) {
636      G4cerr << "G4VUserPhysicsList::PreparePhysicsTable  : No Process Manager for " 
637             << particle->GetParticleName() <<G4endl;
638      G4cerr << particle->GetParticleName() << " should be created in your PhysicsList" <<G4endl;
639      G4Exception("G4VUserPhysicsList::PreparePhysicsTable","No process manager",
640                    FatalException,  particle->GetParticleName() );
641    }
642
643    G4ProcessVector* pVector = pManager->GetProcessList();
644    if (!pVector) {
645      G4cerr << "G4VUserPhysicsList::PreparePhysicsTable  : No Process Vector for " 
646             << particle->GetParticleName() <<G4endl;
647      G4cerr << particle->GetParticleName() << " should be created in your PhysicsList" <<G4endl;
648      G4Exception("G4VUserPhysicsList::PreparePhysicsTable","No process Vector",
649                    FatalException,  particle->GetParticleName() );
650    }
651    for (G4int j=0; j < pVector->size(); ++j) {
652      (*pVector)[j]->PreparePhysicsTable(*particle);
653    }
654  }
655}
656
657///////////////////////////////////////////////////////////////
658void  G4VUserPhysicsList::BuildIntegralPhysicsTable(G4VProcess* process,
659                                                    G4ParticleDefinition* particle)
660{
661  //*********************************************************************
662  // temporary addition to make the integral schema of electromagnetic
663  // processes work.
664  //
665
666  if ( (process->GetProcessName() == "Imsc") ||
667       (process->GetProcessName() == "IeIoni") ||
668       (process->GetProcessName() == "IeBrems") ||
669       (process->GetProcessName() == "Iannihil") ||
670       (process->GetProcessName() == "IhIoni") ||
671       (process->GetProcessName() == "IMuIoni") ||
672       (process->GetProcessName() == "IMuBrems") ||
673       (process->GetProcessName() == "IMuPairProd")  ) {
674#ifdef G4VERBOSE
675    if (verboseLevel>2){
676      G4cout << "G4VUserPhysicsList::BuildIntegralPhysicsTable  ";
677      G4cout << " BuildPhysicsTable is invoked for ";
678      G4cout << process->GetProcessName();
679      G4cout << "(" << particle->GetParticleName() << ")" << G4endl;
680    }
681#endif
682    process->BuildPhysicsTable(*particle);
683  }
684}
685
686///////////////////////////////////////////////////////////////
687void G4VUserPhysicsList::DumpList() const
688{
689  theParticleIterator->reset();
690  G4int idx = 0;
691  while( (*theParticleIterator)() ){
692    G4ParticleDefinition* particle = theParticleIterator->value();
693    G4cout << particle->GetParticleName();
694    if ((idx++ % 4) == 3) {
695      G4cout << G4endl;
696    } else {
697      G4cout << ", ";
698    }     
699  }
700  G4cout << G4endl;
701}
702
703
704
705///////////////////////////////////////////////////////////////
706void G4VUserPhysicsList::DumpCutValuesTable(G4int nParticles) 
707{ 
708  fDisplayThreshold = nParticles; 
709}
710
711///////////////////////////////////////////////////////////////
712void G4VUserPhysicsList::DumpCutValuesTableIfRequested()
713{
714  if(fDisplayThreshold==0) return;
715  G4ProductionCutsTable::GetProductionCutsTable()->DumpCouples();
716  fDisplayThreshold = 0;
717}
718
719
720///////////////////////////////////////////////////////////////
721///////////////////////////////////////////////////////////////
722G4bool G4VUserPhysicsList::StorePhysicsTable(const G4String& directory)
723{
724  G4bool   ascii = fStoredInAscii;
725  G4String dir   = directory;
726  if (dir.isNull()) dir = directoryPhysicsTable; 
727  else directoryPhysicsTable = dir; 
728
729  // store CutsTable info
730  if (!fCutsTable->StoreCutsTable(dir, ascii)) {
731    G4Exception("G4VUserPhysicsList::StorePhysicsTable","Faile to store ",
732                JustWarning,"Cut Table can not be stored");     
733    return false;
734  }
735#ifdef G4VERBOSE 
736  if (verboseLevel>2){
737    G4cout << "G4VUserPhysicsList::StorePhysicsTable   ";
738    G4cout << " Store material and cut values successfully" << G4endl;
739  }
740#endif
741
742  G4bool success= true;
743
744  // loop over all particles in G4ParticleTable
745  theParticleIterator->reset();
746  while( (*theParticleIterator)() ){
747    G4ParticleDefinition* particle = theParticleIterator->value();
748    // Store physics tables for every process for this particle type
749    G4ProcessVector* pVector = (particle->GetProcessManager())->GetProcessList();
750    G4int  j;
751    for ( j=0; j < pVector->size(); ++j) {
752      if (!(*pVector)[j]->StorePhysicsTable(particle,dir,ascii)){   
753        G4String comment =  "Fail to store for " + (*pVector)[j]->GetProcessName();
754        comment += "(" + particle->GetParticleName()  + ")";
755        G4Exception("G4VUserPhysicsList::StorePhysicsTable","Faile to store ",
756                    JustWarning,comment);       
757        success = false;
758      }
759    }
760    // end loop over processes
761  }
762  // end loop over particles
763  return success;
764}
765
766
767
768///////////////////////////////////////////////////////////////
769void  G4VUserPhysicsList::SetPhysicsTableRetrieved(const G4String& directory)
770{
771  fRetrievePhysicsTable = true;
772  if(!directory.isNull()) {
773    directoryPhysicsTable = directory;
774  }
775  fIsCheckedForRetrievePhysicsTable=false;
776  fIsRestoredCutValues = false;
777}
778
779///////////////////////////////////////////////////////////////
780void G4VUserPhysicsList::RetrievePhysicsTable(G4ParticleDefinition* particle, 
781                                              const G4String& directory,
782                                              G4bool          ascii)
783{
784  G4int  j;
785  G4bool success[100]; 
786  // Retrieve physics tables for every process for this particle type
787  G4ProcessVector* pVector = (particle->GetProcessManager())->GetProcessList();
788  for ( j=0; j < pVector->size(); ++j) {
789    success[j] = 
790       (*pVector)[j]->RetrievePhysicsTable(particle,directory,ascii);
791
792    if (!success[j]) {
793#ifdef G4VERBOSE 
794      if (verboseLevel>2){
795        G4cout << "G4VUserPhysicsList::RetrievePhysicsTable   ";
796        G4cout << " Fail to retrieve Physics Table for ";
797        G4cout << (*pVector)[j]->GetProcessName() << G4endl;
798        G4cout << "Calculate Physics Table for " << particle->GetParticleName() << G4endl;
799      }
800#endif
801      (*pVector)[j]->BuildPhysicsTable(*particle);
802    }
803  }
804  for ( j=0; j < pVector->size(); ++j) {
805    // temporary addition to make the integral schema
806    if (!success[j]) BuildIntegralPhysicsTable((*pVector)[j], particle); 
807  }
808}
809
810void G4VUserPhysicsList::ResetCuts()
811{
812#ifdef G4VERBOSE 
813  if (verboseLevel>0){
814    G4cout << "G4VUserPhysicsList::ResetCuts() is obsolete.";
815    G4cout << " This method gives no effect and you can remove it. "<< G4endl;
816  }
817#endif
818}
819
820void G4VUserPhysicsList::SetApplyCuts(G4bool value, const G4String& name)
821{
822#ifdef G4VERBOSE 
823  if (verboseLevel>2){
824    G4cout << "G4VUserPhysicsList::SetApplyCuts for " << name << G4endl;
825  }
826#endif
827  if(name=="all") {
828    theParticleTable->FindParticle("gamma")->SetApplyCutsFlag(value);
829    theParticleTable->FindParticle("e-")->SetApplyCutsFlag(value);
830    theParticleTable->FindParticle("e+")->SetApplyCutsFlag(value);
831    theParticleTable->FindParticle("proton")->SetApplyCutsFlag(value);
832  } else {
833    theParticleTable->FindParticle(name)->SetApplyCutsFlag(value);
834  }
835}
836
837G4bool G4VUserPhysicsList::GetApplyCuts(const G4String& name) const
838{
839  return theParticleTable->FindParticle(name)->GetApplyCutsFlag();
840}
841
842
843
844
845/// obsolete methods
846
847
848void G4VUserPhysicsList::DumpCutValues(const G4String &particle_name)
849{
850  G4cerr << "WARNING !" << G4endl;
851  G4cerr << " Obsolete DumpCutValues() method is invoked for " << particle_name << G4endl;
852  G4cerr << " Please use DumpCutValuesTable() instead." << G4endl;
853  G4cerr << " This dummy method implementation will be removed soon." << G4endl;
854  DumpCutValuesTable();
855}
856
857void G4VUserPhysicsList::DumpCutValues(G4ParticleDefinition* )
858{
859  G4cerr << "WARNING !" << G4endl;
860  G4cerr << " DumpCutValues() became obsolete." << G4endl;
861  G4cerr << " Please use DumpCutValuesTable() instead." << G4endl;
862  G4cerr << " This dummy method implementation will be removed soon." << G4endl;
863  DumpCutValuesTable();
864}
Note: See TracBrowser for help on using the repository browser.