Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4AdjointPhysicsList.cc
이 파일의 문서화 페이지로 가기
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 //
28 //
29 // $Id: G4AdjointPhysicsList.cc 107258 2017-11-07 09:58:16Z gcosmo $
30 //
32 // Class Name: G4AdjointPhysicsList
33 // Author: L. Desorgher
34 // Organisation: SpaceIT GmbH
35 // Contract: ESA contract 21435/08/NL/AT
36 // Customer: ESA/ESTEC
38 
39 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
40 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
41 
42 #include "G4AdjointPhysicsList.hh"
43 #include "G4ProcessManager.hh"
44 #include "G4ParticleTypes.hh"
46 #include "G4SystemOfUnits.hh"
47 
48 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
49 
52  fEminusIonisation(0),fPIonisation(0),
53  fUse_forced_interaction(true),
54  fUse_eionisation(true),fUse_pionisation(true),
55  fUse_brem(true),fUse_compton(true),fUse_ms(true),
56  fUse_egain_fluctuation(true),fUse_peeffect(true),
57  fEmin_adj_models(1.*keV), fEmax_adj_models(1.*MeV),
58  fCS_biasing_factor_compton(1.),fCS_biasing_factor_brem(1.),
59  fCS_biasing_factor_ionisation(1.),fCS_biasing_factor_PEeffect(1.)
60 {
61  defaultCutValue = 1.0*mm;
62  SetVerboseLevel(1);
64 }
65 
66 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
67 
69 {
70 }
72 {
73  // In this method, static member functions should be called
74  // for all particles which you want to use.
75  // This ensures that objects of these particle types will be
76  // created in the program.
82 }
83 
84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
85 
87 {
89 }
90 
91 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
92 
94 {
95  // pseudo-particles
98 
99  // gamma
101 
102  // optical photon
104 }
105 
106 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
107 
109 {
110  // leptons
115 
120 }
121 
122 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
123 
125 {
126 // mesons
138 }
139 
140 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
141 
143 {
144 // barions
149 }
150 
151 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
152 
153 #include"G4AdjointGamma.hh"
154 #include"G4AdjointElectron.hh"
155 #include"G4AdjointProton.hh"
157 {
158 // adjoint_gammma
160 
161 // adjoint_electron
163 
164 // adjoint_proton
166 }
167 
168 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
169 
171 {
173  ConstructEM();
175 }
176 
177 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
178 
179 //#include "G4PEEffectFluoModel.hh"
180 #include "G4ComptonScattering.hh"
181 #include "G4GammaConversion.hh"
182 #include "G4PhotoElectricEffect.hh"
183 #include "G4eMultipleScattering.hh"
185 #include "G4hMultipleScattering.hh"
186 #include "G4eIonisation.hh"
187 #include "G4eBremsstrahlung.hh"
188 #include "G4eplusAnnihilation.hh"
189 #include "G4hIonisation.hh"
190 #include "G4ionIonisation.hh"
191 //#include "G4IonParametrisedLossModel.hh"
192 
193 #include "G4eBremsstrahlung.hh"
195 #include "G4eInverseIonisation.hh"
197 #include "G4AdjointCSManager.hh"
200 #include "G4AdjointComptonModel.hh"
201 #include "G4eInverseCompton.hh"
202 #include "G4InversePEEffect.hh"
205 #include "G4hInverseIonisation.hh"
208 #include "G4IonInverseIonisation.hh"
210 
211 #include "G4AdjointSimManager.hh"
213 
214 #include "G4SystemOfUnits.hh"
215 #include "G4PhysicalConstants.hh"
216 #include "G4UrbanMscModel.hh"
217 #include "G4UrbanAdjointMscModel.hh"
218 #include "G4UrbanAdjointMscModel.hh"
220 
221 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
222 
224 {
225  G4AdjointCSManager* theCSManager =
227  G4AdjointSimManager* theAdjointSimManager =
229 
230  theCSManager->RegisterAdjointParticle(
232 
234  theCSManager->RegisterAdjointParticle(
236 
237  if (fUse_eionisation) {
239  new G4eIonisation();
241  }
242  if (fUse_pionisation) {
245  theCSManager->RegisterAdjointParticle(
247  }
248 
249  G4eBremsstrahlung* theeminusBremsstrahlung = 0;
251  theeminusBremsstrahlung = new G4eBremsstrahlung();
252 
253  G4ComptonScattering* theComptonScattering =0;
254  if (fUse_compton) theComptonScattering = new G4ComptonScattering();
255 
256  G4PhotoElectricEffect* thePEEffect =0;
257  if (fUse_peeffect) thePEEffect = new G4PhotoElectricEffect();
258 
259  G4eMultipleScattering* theeminusMS = 0;
260  G4hMultipleScattering* thepMS= 0;
261  G4eAdjointMultipleScattering* theeminusAdjointMS = 0;
262  if (fUse_ms) {
263  theeminusMS = new G4eMultipleScattering();
264  G4UrbanMscModel* msc1 = new G4UrbanMscModel();
265  theeminusMS->SetEmModel(msc1);
266  theeminusAdjointMS = new G4eAdjointMultipleScattering();
268  theeminusAdjointMS->SetEmModel(msc2);
269  thepMS = new G4hMultipleScattering();
270  }
271 
272  G4VProcess* theGammaConversion =0;
273  if (fUse_gamma_conversion) theGammaConversion = new G4GammaConversion();
274  //Define adjoint e- ionisation
275  //-------------------
276  G4AdjointeIonisationModel* theeInverseIonisationModel = 0;
277  G4eInverseIonisation* theeInverseIonisationProjToProjCase = 0 ;
278  G4eInverseIonisation* theeInverseIonisationProdToProjCase = 0;
279  if (fUse_eionisation) {
280  theeInverseIonisationModel = new G4AdjointeIonisationModel();
281  theeInverseIonisationModel->SetHighEnergyLimit(
283  theeInverseIonisationModel->SetLowEnergyLimit(
285  theeInverseIonisationModel->SetCSBiasingFactor(
287  theeInverseIonisationProjToProjCase =
288  new G4eInverseIonisation(true,"Inv_eIon",
289  theeInverseIonisationModel);
290  theeInverseIonisationProdToProjCase =
291  new G4eInverseIonisation(false,"Inv_eIon1",
292  theeInverseIonisationModel);
293  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("e-"));
294  }
295 
296  //Define adjoint Bremsstrahlung
297  //-------------------------------
298  G4AdjointBremsstrahlungModel* theeInverseBremsstrahlungModel = 0;
299  G4eInverseBremsstrahlung* theeInverseBremsstrahlungProjToProjCase = 0;
300  G4eInverseBremsstrahlung* theeInverseBremsstrahlungProdToProjCase = 0;
301  G4AdjointForcedInteractionForGamma* theForcedInteractionForGamma = 0;
302  if (fUse_brem && fUse_eionisation) {
303  theeInverseBremsstrahlungModel = new G4AdjointBremsstrahlungModel();
304  theeInverseBremsstrahlungModel->SetHighEnergyLimit(fEmax_adj_models*1.01);
305  theeInverseBremsstrahlungModel->SetLowEnergyLimit(fEmin_adj_models);
306  theeInverseBremsstrahlungModel->SetCSBiasingFactor(
308  theeInverseBremsstrahlungProjToProjCase = new G4eInverseBremsstrahlung(
309  true,"Inv_eBrem",theeInverseBremsstrahlungModel);
310  theeInverseBremsstrahlungProdToProjCase = new G4eInverseBremsstrahlung(
311  false,"Inv_eBrem1",theeInverseBremsstrahlungModel);
312  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("e-"));
313  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("gamma"));
314 
315  if (!fUse_forced_interaction) theeInverseBremsstrahlungProdToProjCase
316  = new G4eInverseBremsstrahlung(false,G4String("Inv_eBrem1"),
317  theeInverseBremsstrahlungModel);
318  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("e-"));
319  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("gamma"));
321  theForcedInteractionForGamma =
322  new G4AdjointForcedInteractionForGamma("ReverseGammaForcedInteraction");
323  theForcedInteractionForGamma->RegisterAdjointBremModel(
324  theeInverseBremsstrahlungModel);
325  }
326  }
327 
328 
329  //Define adjoint Compton
330  //---------------------
331 
332  G4AdjointComptonModel* theeInverseComptonModel = 0;
333  G4eInverseCompton* theeInverseComptonProjToProjCase = 0;
334  G4eInverseCompton* theeInverseComptonProdToProjCase = 0;
335 
336  if (fUse_compton) {
337  theeInverseComptonModel = new G4AdjointComptonModel();
338  theeInverseComptonModel->SetHighEnergyLimit(fEmax_adj_models);
339  theeInverseComptonModel->SetLowEnergyLimit(fEmin_adj_models);
340  theeInverseComptonModel->SetDirectProcess(theComptonScattering);
341  theeInverseComptonModel->SetUseMatrix(false);
342 
343  theeInverseComptonModel->SetCSBiasingFactor( fCS_biasing_factor_compton);
344  if (!fUse_forced_interaction) theeInverseComptonProjToProjCase =
345  new G4eInverseCompton(true,"Inv_Compt",theeInverseComptonModel);
346  theeInverseComptonProdToProjCase = new G4eInverseCompton(false,"Inv_Compt1",
347  theeInverseComptonModel);
349  if (!theForcedInteractionForGamma ) theForcedInteractionForGamma =
350  new G4AdjointForcedInteractionForGamma("ReverseGammaForcedInteraction");
351  theForcedInteractionForGamma->
352  RegisterAdjointComptonModel(theeInverseComptonModel);
353  }
354  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("e-"));
355  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("gamma"));
356  }
357 
358  //Define adjoint PEEffect
359  //---------------------
360  G4AdjointPhotoElectricModel* theInversePhotoElectricModel = 0;
361  G4InversePEEffect* theInversePhotoElectricProcess = 0;
362 
363  if (fUse_peeffect) {
364  theInversePhotoElectricModel = new G4AdjointPhotoElectricModel();
365  theInversePhotoElectricModel->SetHighEnergyLimit(fEmax_adj_models);
366  theInversePhotoElectricModel->SetLowEnergyLimit(fEmin_adj_models);
367  theInversePhotoElectricModel->SetCSBiasingFactor(
369  theInversePhotoElectricProcess = new G4InversePEEffect("Inv_PEEffect",
370  theInversePhotoElectricModel);
371  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("e-"));
372  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("gamma"));
373  }
374 
375 
376  //Define adjoint ionisation for protons
377  //---------------------
378  G4AdjointhIonisationModel* thepInverseIonisationModel = 0;
379  G4hInverseIonisation* thepInverseIonisationProjToProjCase = 0 ;
380  G4hInverseIonisation* thepInverseIonisationProdToProjCase = 0;
381  if (fUse_pionisation) {
382  thepInverseIonisationModel = new G4AdjointhIonisationModel(
383  G4Proton::Proton());
384  thepInverseIonisationModel->SetHighEnergyLimit(fEmax_adj_models);
385  thepInverseIonisationModel->SetLowEnergyLimit(fEmin_adj_models);
386  thepInverseIonisationModel->SetUseMatrix(false);
387  thepInverseIonisationProjToProjCase = new G4hInverseIonisation(true,
388  "Inv_pIon",thepInverseIonisationModel);
389  thepInverseIonisationProdToProjCase = new G4hInverseIonisation(false,
390  "Inv_pIon1",thepInverseIonisationModel);
391  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("e-"));
392  theAdjointSimManager->ConsiderParticleAsPrimary(G4String("proton"));
393  }
394 
395  //Declare the processes active for the different particles
396  //--------------------------------------------------------
398  particleIterator->reset();
399  while( (*particleIterator)() ){
400  G4ParticleDefinition* particle = particleIterator->value();
401  G4ProcessManager* pmanager = particle->GetProcessManager();
402  if (!pmanager) {
403  pmanager = new G4ProcessManager(particle);
404  particle->SetProcessManager(pmanager);
405  }
406 
407  G4String particleName = particle->GetParticleName();
408  if (particleName == "e-") {
409  if (fUse_ms && fUse_eionisation) pmanager->AddProcess(theeminusMS);
410  if (fUse_eionisation){
411  pmanager->AddProcess(fEminusIonisation);
413  RegisterEnergyLossProcess(fEminusIonisation,particle);
414  }
415  if (fUse_brem && fUse_eionisation) {
416  pmanager->AddProcess(theeminusBremsstrahlung);
418  RegisterEnergyLossProcess(theeminusBremsstrahlung,particle);
419  }
420  G4int n_order=0;
421  if (fUse_ms && fUse_eionisation) {
422  n_order++;
423  pmanager->SetProcessOrdering(theeminusMS, idxAlongStep,n_order);
424  }
425  if (fUse_eionisation) {
426  n_order++;
428  }
429  if (fUse_brem && fUse_eionisation) {
430  n_order++;
431  pmanager->SetProcessOrdering(theeminusBremsstrahlung,
432  idxAlongStep,n_order);
433  }
434  n_order=0;
435  if (fUse_ms && fUse_eionisation) {
436  n_order++;
437  pmanager->SetProcessOrdering(theeminusMS,idxPostStep,n_order);
438  }
439  if (fUse_eionisation) {
440  n_order++;
442  }
443  if (fUse_brem && fUse_eionisation) {
444  n_order++;
445  pmanager->SetProcessOrdering(theeminusBremsstrahlung,idxPostStep,
446  n_order);
447  }
448  }
449 
450  if (particleName == "adj_e-") {
451  G4ContinuousGainOfEnergy* theContinuousGainOfEnergy =0;
452  if (fUse_eionisation ) {
453  theContinuousGainOfEnergy= new G4ContinuousGainOfEnergy();
454  theContinuousGainOfEnergy->SetLossFluctuations(
456  theContinuousGainOfEnergy->SetDirectEnergyLossProcess(
458  theContinuousGainOfEnergy->SetDirectParticle(G4Electron::Electron());
459  pmanager->AddProcess(theContinuousGainOfEnergy);
460  }
461  G4int n_order=0;
462  if (fUse_ms) {
463  n_order++;
464  pmanager->AddProcess(theeminusAdjointMS);
465  pmanager->SetProcessOrdering(theeminusAdjointMS,
466  idxAlongStep,n_order);
467  }
468  n_order++;
469  pmanager->SetProcessOrdering(theContinuousGainOfEnergy,idxAlongStep,
470  n_order);
471 
472  n_order++;
473  G4AdjointAlongStepWeightCorrection* theAlongStepWeightCorrection =
475  pmanager->AddProcess(theAlongStepWeightCorrection);
476  pmanager->SetProcessOrdering(theAlongStepWeightCorrection,
477  idxAlongStep,
478  n_order);
479  n_order=0;
480  if (fUse_eionisation) {
481  pmanager->AddProcess(theeInverseIonisationProjToProjCase);
482  pmanager->AddProcess(theeInverseIonisationProdToProjCase);
483  n_order++;
484  pmanager->SetProcessOrdering(theeInverseIonisationProjToProjCase,
485  idxPostStep,n_order);
486  n_order++;
487  pmanager->SetProcessOrdering(theeInverseIonisationProdToProjCase,
488  idxPostStep,n_order);
489  }
490  if (fUse_brem && fUse_eionisation) {
491  pmanager->AddProcess(theeInverseBremsstrahlungProjToProjCase);
492  n_order++;
493  pmanager->SetProcessOrdering(
494  theeInverseBremsstrahlungProjToProjCase,
495  idxPostStep,n_order);
496  }
497 
498  if (fUse_compton) {
499  pmanager->AddProcess(theeInverseComptonProdToProjCase);
500  n_order++;
501  pmanager->SetProcessOrdering(theeInverseComptonProdToProjCase,
502  idxPostStep,n_order);
503  }
504  if (fUse_peeffect) {
505  pmanager->AddDiscreteProcess(theInversePhotoElectricProcess);
506  n_order++;
507  pmanager->SetProcessOrdering(theInversePhotoElectricProcess,
508  idxPostStep,n_order);
509  }
510  if (fUse_pionisation) {
511  pmanager->AddProcess(thepInverseIonisationProdToProjCase);
512  n_order++;
513  pmanager->SetProcessOrdering(thepInverseIonisationProdToProjCase,
514  idxPostStep,n_order);
515  }
516  if (fUse_ms && fUse_eionisation) {
517  n_order++;
518  pmanager->SetProcessOrdering(theeminusAdjointMS,
519  idxPostStep,n_order);
520  }
521  }
522 
523 
524  if(particleName == "adj_gamma") {
525  G4int n_order=0;
527  G4AdjointAlongStepWeightCorrection* theAlongStepWeightCorrection =
529  pmanager->AddProcess(theAlongStepWeightCorrection);
530  pmanager->SetProcessOrdering(theAlongStepWeightCorrection,
531  idxAlongStep,1);
532 
533  if (fUse_brem && fUse_eionisation) {
534  pmanager->AddProcess(theeInverseBremsstrahlungProdToProjCase);
535  n_order++;
536  pmanager->SetProcessOrdering(
537  theeInverseBremsstrahlungProdToProjCase,
538  idxPostStep,n_order);
539  }
540  if (fUse_compton) {
541  pmanager->AddDiscreteProcess(theeInverseComptonProjToProjCase);
542  n_order++;
543  pmanager->SetProcessOrdering(theeInverseComptonProjToProjCase,
544  idxPostStep,n_order);
545  }
546  }
547  else {
548  if (theForcedInteractionForGamma) {
549  pmanager->AddProcess(theForcedInteractionForGamma);
550  n_order++;
551  pmanager->SetProcessOrdering(theForcedInteractionForGamma,
552  idxPostStep,n_order);
553  pmanager->SetProcessOrdering(theForcedInteractionForGamma,
554  idxAlongStep,n_order);
555  }
556  }
557  }
558 
559  if (particleName == "gamma") {
560  if (fUse_compton) {
561  pmanager->AddDiscreteProcess(theComptonScattering);
563  RegisterEmProcess(theComptonScattering,particle);
564  }
565  if (fUse_peeffect) {
566  pmanager->AddDiscreteProcess(thePEEffect);
568  RegisterEmProcess(thePEEffect,particle);
569  }
570  if (fUse_gamma_conversion) {
571  pmanager->AddDiscreteProcess(theGammaConversion);
572  }
573  }
574 
575  if (particleName == "e+" && fUse_gamma_conversion) {//positron
576  G4VProcess* theeplusMultipleScattering = new G4eMultipleScattering();
577  G4VProcess* theeplusIonisation = new G4eIonisation();
578  G4VProcess* theeplusBremsstrahlung = new G4eBremsstrahlung();
579  G4VProcess* theeplusAnnihilation = new G4eplusAnnihilation();
580 
581  // add processes
582  pmanager->AddProcess(theeplusMultipleScattering);
583  pmanager->AddProcess(theeplusIonisation);
584  pmanager->AddProcess(theeplusBremsstrahlung);
585  pmanager->AddProcess(theeplusAnnihilation);
586 
587  // set ordering for AtRestDoIt
588  pmanager->SetProcessOrderingToFirst(theeplusAnnihilation, idxAtRest);
589 
590  // set ordering for AlongStepDoIt
591  pmanager->SetProcessOrdering(theeplusMultipleScattering,
592  idxAlongStep,1);
593  pmanager->SetProcessOrdering(theeplusIonisation, idxAlongStep,2);
594  pmanager->SetProcessOrdering(theeplusBremsstrahlung,idxAlongStep,3);
595 
596  // set ordering for PostStepDoIt
597  pmanager->SetProcessOrdering(theeplusMultipleScattering,
598  idxPostStep,1);
599  pmanager->SetProcessOrdering(theeplusIonisation,idxPostStep,2);
600  pmanager->SetProcessOrdering(theeplusBremsstrahlung,idxPostStep,3);
601  pmanager->SetProcessOrdering(theeplusAnnihilation,idxPostStep,4);
602  }
603  if (particleName == "proton" && fUse_pionisation) {
604  if (fUse_ms && fUse_pionisation) pmanager->AddProcess(thepMS);
605 
606  if (fUse_pionisation){
607  pmanager->AddProcess(fPIonisation);
609  RegisterEnergyLossProcess(fPIonisation,particle);
610  }
611 
612  G4int n_order=0;
613  if (fUse_ms && fUse_pionisation) {
614  n_order++;
615  pmanager->SetProcessOrdering(thepMS, idxAlongStep,n_order);
616  }
617 
618  if (fUse_pionisation) {
619  n_order++;
620  pmanager->SetProcessOrdering(fPIonisation,idxAlongStep,n_order);
621  }
622 
623  n_order=0;
624  if (fUse_ms && fUse_pionisation) {
625  n_order++;
626  pmanager->SetProcessOrdering(thepMS, idxPostStep,n_order);
627  }
628 
629  if (fUse_pionisation) {
630  n_order++;
631  pmanager->SetProcessOrdering(fPIonisation,idxPostStep,n_order);
632  }
633 
634  }
635 
636  if (particleName == "adj_proton" && fUse_pionisation) {
637  G4ContinuousGainOfEnergy* theContinuousGainOfEnergy =0;
638  if (fUse_pionisation ) {
639  theContinuousGainOfEnergy= new G4ContinuousGainOfEnergy();
640  theContinuousGainOfEnergy->SetLossFluctuations(
642  theContinuousGainOfEnergy->SetDirectEnergyLossProcess(fPIonisation);
643  theContinuousGainOfEnergy->SetDirectParticle(G4Proton::Proton());
644  pmanager->AddProcess(theContinuousGainOfEnergy);
645  }
646 
647  G4int n_order=0;
648  if (fUse_ms) {
649  n_order++;
650  pmanager->AddProcess(thepMS);
651  pmanager->SetProcessOrdering(thepMS, idxAlongStep,n_order);
652  }
653 
654  n_order++;
655  pmanager->SetProcessOrdering(theContinuousGainOfEnergy,idxAlongStep,
656  n_order);
657 
658  n_order++;
659  G4AdjointAlongStepWeightCorrection* theAlongStepWeightCorrection =
661  pmanager->AddProcess(theAlongStepWeightCorrection);
662  pmanager->SetProcessOrdering(theAlongStepWeightCorrection,
663  idxAlongStep,
664  n_order);
665  n_order=0;
666  if (fUse_pionisation) {
667  pmanager->AddProcess(thepInverseIonisationProjToProjCase);
668  n_order++;
669  pmanager->SetProcessOrdering(
670  thepInverseIonisationProjToProjCase,
671  idxPostStep,n_order);
672  }
673 
674  if (fUse_ms && fUse_pionisation) {
675  n_order++;
676  pmanager->SetProcessOrdering(thepMS,idxPostStep,n_order);
677  }
678  }
679  }
680 }
681 
682 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
683 
684 #include "G4Decay.hh"
686 {
687  // Add Decay Process
688  G4Decay* theDecayProcess = new G4Decay();
690  particleIterator->reset();
691  while( (*particleIterator)() ){
692  G4ParticleDefinition* particle = particleIterator->value();
693  G4ProcessManager* pmanager = particle->GetProcessManager();
694  if (theDecayProcess->IsApplicable(*particle)) {
695  pmanager ->AddProcess(theDecayProcess);
696  // set ordering for PostStepDoIt and AtRestDoIt
697  pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
698  pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
699  }
700  }
701 }
702 
703 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
704 
706 {
707  if (verboseLevel >0){
708  G4cout << "G4AdjointPhysicsList::SetCuts:";
709  G4cout << "CutLength : " << G4BestUnit(defaultCutValue,"Length") << G4endl;
710  }
711 
712  // set cut values for gamma at first and for e- second and next for e+,
713  // because some processes for e+/e- need cut values for gamma
714  //
715  SetCutValue(defaultCutValue, "gamma");
718 
720 }
721 
722 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
Definition of the G4AdjointPhysicsMessenger class.
static constexpr double MeV
Definition: G4SIunits.hh:214
static G4AntiNeutron * AntiNeutronDefinition()
void SetDirectParticle(G4ParticleDefinition *p)
G4int AddDiscreteProcess(G4VProcess *aProcess, G4int ord=ordDefault)
static G4ChargedGeantino * ChargedGeantinoDefinition()
static constexpr double keV
Definition: G4SIunits.hh:216
static constexpr double mm
Definition: G4SIunits.hh:115
static G4KaonZero * KaonZeroDefinition()
Definition: G4KaonZero.cc:99
#define G4endl
Definition: G4ios.hh:61
static G4AdjointElectron * AdjointElectronDefinition()
const G4String & GetParticleName() const
static G4AdjointElectron * AdjointElectron()
virtual G4bool IsApplicable(const G4ParticleDefinition &) override
Definition: G4Decay.cc:89
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4PionPlus * PionPlusDefinition()
Definition: G4PionPlus.cc:93
static G4KaonZeroShort * KaonZeroShortDefinition()
G4eIonisation * fEminusIonisation
void SetProcessOrdering(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt, G4int ordDoIt=ordDefault)
static G4AdjointProton * AdjointProtonDefinition()
static G4MuonMinus * MuonMinusDefinition()
Definition: G4MuonMinus.cc:95
static G4AdjointGamma * AdjointGamma()
G4ParticleTable::G4PTblDicIterator * GetParticleIterator() const
G4int AddProcess(G4VProcess *aProcess, G4int ordAtRestDoIt=ordInActive, G4int ordAlongSteptDoIt=ordInActive, G4int ordPostStepDoIt=ordInActive)
void RegisterAdjointBremModel(G4VEmAdjointModel *aAdjointBremModel)
void SetProcessManager(G4ProcessManager *aProcessManager)
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
static G4AdjointGamma * AdjointGammaDefinition()
static G4PionZero * PionZeroDefinition()
Definition: G4PionZero.cc:103
void SetProcessOrderingToFirst(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt)
static G4OpticalPhoton * OpticalPhotonDefinition()
static G4AdjointProton * AdjointProton()
void RegisterAdjointParticle(G4ParticleDefinition *aPartDef)
void SetDirectEnergyLossProcess(G4VEnergyLossProcess *aProcess)
static G4PionMinus * PionMinusDefinition()
Definition: G4PionMinus.cc:93
static G4Electron * Electron()
Definition: G4Electron.cc:94
static G4Electron * ElectronDefinition()
Definition: G4Electron.cc:89
static G4AdjointCSManager * GetAdjointCSManager()
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
static G4AntiNeutrinoMu * AntiNeutrinoMuDefinition()
static G4KaonMinus * KaonMinusDefinition()
Definition: G4KaonMinus.cc:108
static G4Gamma * GammaDefinition()
Definition: G4Gamma.cc:81
int G4int
Definition: G4Types.hh:78
void SetEmModel(G4VMscModel *, size_t index=0)
G4ProcessManager * GetProcessManager() const
static G4NeutrinoE * NeutrinoEDefinition()
Definition: G4NeutrinoE.cc:80
void SetLossFluctuations(G4bool val)
Definition of the G4AdjointPhysicsList class.
virtual void ConstructParticle()
static G4KaonZeroLong * KaonZeroLongDefinition()
virtual void SetCSBiasingFactor(G4double aVal)
static G4Neutron * NeutronDefinition()
Definition: G4Neutron.cc:99
std::vector< G4InuclElementaryParticle >::iterator particleIterator
Definition: G4BigBanger.cc:65
void SetLowEnergyLimit(G4double aVal)
G4GLOB_DLL std::ostream G4cout
static G4AntiKaonZero * AntiKaonZeroDefinition()
static G4Eta * EtaDefinition()
Definition: G4Eta.cc:104
void SetCutValue(G4double aCut, const G4String &pname)
G4AdjointPhysicsMessenger * fPhysicsMessenger
static G4Positron * PositronDefinition()
Definition: G4Positron.cc:89
static G4NeutrinoMu * NeutrinoMuDefinition()
Definition: G4NeutrinoMu.cc:80
static G4Geantino * GeantinoDefinition()
Definition: G4Geantino.cc:82
void SetHighEnergyLimit(G4double aVal)
static G4KaonPlus * KaonPlusDefinition()
Definition: G4KaonPlus.cc:108
static G4MuonPlus * MuonPlusDefinition()
Definition: G4MuonPlus.cc:94
void SetVerboseLevel(G4int value)
static G4AdjointSimManager * GetInstance()
void ConsiderParticleAsPrimary(const G4String &particle_name)
void DumpCutValuesTable(G4int flag=1)
static G4AntiProton * AntiProtonDefinition()
Definition: G4AntiProton.cc:88
void SetDirectProcess(G4VEmProcess *aProcess)
static G4AntiNeutrinoE * AntiNeutrinoEDefinition()
static G4EtaPrime * EtaPrimeDefinition()
Definition: G4EtaPrime.cc:106
G4hIonisation * fPIonisation
void SetLossFluctuationFlag(bool aBool)
void SetUseMatrix(G4bool aBool)