00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055 #include "G4EmProcessOptions.hh"
00056 #include "G4SystemOfUnits.hh"
00057 #include "G4LossTableManager.hh"
00058 #include "G4VEmProcess.hh"
00059 #include "G4VEnergyLossProcess.hh"
00060 #include "G4VMultipleScattering.hh"
00061 #include "G4Region.hh"
00062 #include "G4RegionStore.hh"
00063 #include "G4VAtomDeexcitation.hh"
00064
00065
00066
00067 G4EmProcessOptions::G4EmProcessOptions()
00068 {
00069 theManager = G4LossTableManager::Instance();
00070 }
00071
00072
00073
00074 G4EmProcessOptions::~G4EmProcessOptions()
00075 {
00076 }
00077
00078
00079
00080 void G4EmProcessOptions::SetLossFluctuations(G4bool val)
00081 {
00082 theManager->SetLossFluctuations(val);
00083 }
00084
00085
00086
00087 void G4EmProcessOptions::SetSubCutoff(G4bool val, const G4Region* r)
00088 {
00089 theManager->SetSubCutoff(val, r);
00090 }
00091
00092
00093
00094 void G4EmProcessOptions::SetIntegral(G4bool val)
00095 {
00096 theManager->SetIntegral(val);
00097 }
00098
00099
00100
00101 void G4EmProcessOptions::SetMinSubRange(G4double val)
00102 {
00103 theManager->SetMinSubRange(val);
00104 }
00105
00106
00107
00108 void G4EmProcessOptions::SetMinEnergy(G4double val)
00109 {
00110 theManager->SetMinEnergy(val);
00111 }
00112
00113
00114
00115 void G4EmProcessOptions::SetMaxEnergy(G4double val)
00116 {
00117 theManager->SetMaxEnergy(val);
00118 }
00119
00120
00121
00122 void G4EmProcessOptions::SetMaxEnergyForCSDARange(G4double val)
00123 {
00124 theManager->SetMaxEnergyForCSDARange(val);
00125 }
00126
00127
00128
00129 void G4EmProcessOptions::SetMaxEnergyForMuons(G4double val)
00130 {
00131 theManager->SetMaxEnergyForMuons(val);
00132 }
00133
00134
00135
00136 void G4EmProcessOptions::SetDEDXBinning(G4int val)
00137 {
00138 theManager->SetDEDXBinning(val);
00139 }
00140
00141
00142
00143 void G4EmProcessOptions::SetDEDXBinningForCSDARange(G4int val)
00144 {
00145 theManager->SetDEDXBinningForCSDARange(val);
00146 }
00147
00148
00149
00150 void G4EmProcessOptions::SetLambdaBinning(G4int val)
00151 {
00152 theManager->SetLambdaBinning(val);
00153 }
00154
00155
00156
00157 void G4EmProcessOptions::SetStepFunction(G4double v1, G4double v2)
00158 {
00159 theManager->SetStepFunction(v1, v2);
00160 }
00161
00162
00163
00164 void G4EmProcessOptions::SetRandomStep(G4bool val)
00165 {
00166 theManager->SetRandomStep(val);
00167 }
00168
00169
00170
00171 void G4EmProcessOptions::SetApplyCuts(G4bool val)
00172 {
00173 const std::vector<G4VEmProcess*>& w =
00174 theManager->GetEmProcessVector();
00175 std::vector<G4VEmProcess*>::const_iterator itp;
00176 for(itp = w.begin(); itp != w.end(); itp++) {
00177 G4VEmProcess* q = *itp;
00178 if(q) { q->SetApplyCuts(val); }
00179 }
00180 }
00181
00182
00183
00184 void G4EmProcessOptions::SetBuildCSDARange(G4bool val)
00185 {
00186 theManager->SetBuildCSDARange(val);
00187 }
00188
00189
00190
00191 void G4EmProcessOptions::SetVerbose(G4int val, const G4String& name)
00192 {
00193 G4bool all = false;
00194 if("all" == name) { all = true; }
00195 const std::vector<G4VEnergyLossProcess*>& v =
00196 theManager->GetEnergyLossProcessVector();
00197
00198 if(all) {
00199 theManager->SetVerbose(val);
00200 return;
00201 }
00202
00203 std::vector<G4VEnergyLossProcess*>::const_iterator itr;
00204 for(itr = v.begin(); itr != v.end(); ++itr) {
00205 G4VEnergyLossProcess* p = *itr;
00206 if(p) {
00207 if (p->GetProcessName() == name) { p->SetVerboseLevel(val); }
00208 }
00209 }
00210 const std::vector<G4VEmProcess*>& w =
00211 theManager->GetEmProcessVector();
00212 std::vector<G4VEmProcess*>::const_iterator itp;
00213 for(itp = w.begin(); itp != w.end(); itp++) {
00214 G4VEmProcess* q = *itp;
00215 if(q) {
00216 if (q->GetProcessName() == name) { q->SetVerboseLevel(val); }
00217 }
00218 }
00219 const std::vector<G4VMultipleScattering*>& u =
00220 theManager->GetMultipleScatteringVector();
00221 std::vector<G4VMultipleScattering*>::const_iterator itm;
00222 for(itm = u.begin(); itm != u.end(); itm++) {
00223 G4VMultipleScattering* msc = *itm;
00224 if(s) {
00225 if (msc->GetProcessName() == name) { msc->SetVerboseLevel(val); }
00226 }
00227 }
00228 }
00229
00230
00231
00232 void G4EmProcessOptions::SetLambdaFactor(G4double val)
00233 {
00234 const std::vector<G4VEnergyLossProcess*>& v =
00235 theManager->GetEnergyLossProcessVector();
00236 std::vector<G4VEnergyLossProcess*>::const_iterator itr;
00237 for(itr = v.begin(); itr != v.end(); itr++) {
00238 G4VEnergyLossProcess* p = *itr;
00239 if(p) { p->SetLambdaFactor(val); }
00240 }
00241 }
00242
00243
00244
00245 void
00246 G4EmProcessOptions::SetDeexcitationActiveRegion(const G4String& rname,
00247 G4bool valDeexcitation,
00248 G4bool valAuger,
00249 G4bool valPIXE)
00250 {
00251 G4VAtomDeexcitation* ad = theManager-> AtomDeexcitation();
00252 if(ad) {
00253 ad->SetDeexcitationActiveRegion(rname, valDeexcitation,
00254 valAuger,valPIXE);
00255 }
00256 }
00257
00258
00259
00260 void G4EmProcessOptions::SetFluo(G4bool val)
00261 {
00262 G4VAtomDeexcitation* ad = theManager-> AtomDeexcitation();
00263 if(ad) { ad->SetFluo(val); }
00264 }
00265
00266
00267
00268 void G4EmProcessOptions::SetAuger(G4bool val)
00269 {
00270 G4VAtomDeexcitation* ad = theManager-> AtomDeexcitation();
00271 if(ad) { ad->SetAuger(val); }
00272 }
00273
00274
00275
00276 void G4EmProcessOptions::SetPIXE(G4bool val)
00277 {
00278 G4VAtomDeexcitation* ad = theManager-> AtomDeexcitation();
00279 if(ad) { ad->SetPIXE(val); }
00280 }
00281
00282
00283
00284 void G4EmProcessOptions::SetPIXECrossSectionModel(const G4String& mname)
00285 {
00286 G4VAtomDeexcitation* ad = theManager-> AtomDeexcitation();
00287 if(ad) { ad->SetPIXECrossSectionModel(mname); }
00288 }
00289
00290
00291
00292 void
00293 G4EmProcessOptions::SetPIXEElectronCrossSectionModel(const G4String& mname)
00294 {
00295 G4VAtomDeexcitation* ad = theManager-> AtomDeexcitation();
00296 if(ad) { ad->SetPIXEElectronCrossSectionModel(mname); }
00297 }
00298
00299
00300
00301 void G4EmProcessOptions::SetMscStepLimitation(G4MscStepLimitType val)
00302 {
00303 const std::vector<G4VMultipleScattering*>& u =
00304 theManager->GetMultipleScatteringVector();
00305 std::vector<G4VMultipleScattering*>::const_iterator itm;
00306 for(itm = u.begin(); itm != u.end(); itm++) {
00307 if(*itm) (*itm)->SetStepLimitType(val);
00308 }
00309 }
00310
00311
00312
00313 void G4EmProcessOptions::SetMscLateralDisplacement(G4bool val)
00314 {
00315 const std::vector<G4VMultipleScattering*>& u =
00316 theManager->GetMultipleScatteringVector();
00317 std::vector<G4VMultipleScattering*>::const_iterator itm;
00318 for(itm = u.begin(); itm != u.end(); itm++) {
00319 if(*itm) { (*itm)->SetLateralDisplasmentFlag(val); }
00320 }
00321 }
00322
00323
00324
00325 void G4EmProcessOptions::SetSkin(G4double val)
00326 {
00327 if(val < 0.0) return;
00328 const std::vector<G4VMultipleScattering*>& u =
00329 theManager->GetMultipleScatteringVector();
00330 std::vector<G4VMultipleScattering*>::const_iterator itm;
00331 for(itm = u.begin(); itm != u.end(); itm++) {
00332 if(*itm) { (*itm)->SetSkin(val); }
00333 }
00334 }
00335
00336
00337
00338 void G4EmProcessOptions::SetMscRangeFactor(G4double val)
00339 {
00340 if(val < 0.0) return;
00341 const std::vector<G4VMultipleScattering*>& u =
00342 theManager->GetMultipleScatteringVector();
00343 std::vector<G4VMultipleScattering*>::const_iterator itm;
00344 for(itm = u.begin(); itm != u.end(); itm++) {
00345 if(*itm) { (*itm)->SetRangeFactor(val); }
00346 }
00347 }
00348
00349
00350
00351 void G4EmProcessOptions::SetMscGeomFactor(G4double val)
00352 {
00353 if(val < 0.0) { return; }
00354 const std::vector<G4VMultipleScattering*>& u =
00355 theManager->GetMultipleScatteringVector();
00356 std::vector<G4VMultipleScattering*>::const_iterator itm;
00357 for(itm = u.begin(); itm != u.end(); itm++) {
00358 if(*itm) { (*itm)->SetGeomFactor(val); }
00359 }
00360 }
00361
00362
00363
00364 void G4EmProcessOptions::SetPolarAngleLimit(G4double val)
00365 {
00366 const std::vector<G4VMultipleScattering*>& u =
00367 theManager->GetMultipleScatteringVector();
00368 std::vector<G4VMultipleScattering*>::const_iterator itm;
00369 for(itm = u.begin(); itm != u.end(); itm++) {
00370 if(*itm) { (*itm)->SetPolarAngleLimit(val); }
00371 }
00372 const std::vector<G4VEmProcess*>& w =
00373 theManager->GetEmProcessVector();
00374 std::vector<G4VEmProcess*>::const_iterator itp;
00375 for(itp = w.begin(); itp != w.end(); itp++) {
00376 G4VEmProcess* q = *itp;
00377 if(q) { q->SetPolarAngleLimit(val); }
00378 }
00379 }
00380
00381
00382
00383 void G4EmProcessOptions::SetFactorForAngleLimit(G4double val)
00384 {
00385 theManager->SetFactorForAngleLimit(val);
00386 }
00387
00388
00389
00390 void G4EmProcessOptions::SetLPMFlag(G4bool val)
00391 {
00392 theManager->SetLPMFlag(val);
00393 }
00394
00395
00396
00397 void G4EmProcessOptions::SetSplineFlag(G4bool val)
00398 {
00399 theManager->SetSplineFlag(val);
00400 }
00401
00402
00403
00404 void G4EmProcessOptions::SetLinearLossLimit(G4double val)
00405 {
00406 theManager->SetLinearLossLimit(val);
00407 }
00408
00409
00410
00411 void G4EmProcessOptions::SetBremsstrahlungTh(G4double val)
00412 {
00413 theManager->SetBremsstrahlungTh(val);
00414 }
00415
00416
00417
00418 void
00419 G4EmProcessOptions::SetProcessBiasingFactor(const G4String& name, G4double val,
00420 G4bool flag)
00421 {
00422 const std::vector<G4VEnergyLossProcess*>& v =
00423 theManager->GetEnergyLossProcessVector();
00424 std::vector<G4VEnergyLossProcess*>::const_iterator itr;
00425 for(itr = v.begin(); itr != v.end(); ++itr) {
00426 G4VEnergyLossProcess* p = *itr;
00427 if(p) {
00428 if (p->GetProcessName() == name) {
00429 p->SetCrossSectionBiasingFactor(val, flag);
00430 }
00431 }
00432 }
00433 const std::vector<G4VEmProcess*>& w =
00434 theManager->GetEmProcessVector();
00435 std::vector<G4VEmProcess*>::const_iterator itp;
00436 for(itp = w.begin(); itp != w.end(); itp++) {
00437 G4VEmProcess* q = *itp;
00438 if(q) {
00439 if (q->GetProcessName() == name) {
00440 q->SetCrossSectionBiasingFactor(val, flag);
00441 }
00442 }
00443 }
00444 }
00445
00446
00447
00448 void
00449 G4EmProcessOptions::ActivateForcedInteraction(const G4String& name,
00450 G4double length,
00451 const G4String& region,
00452 G4bool flag)
00453 {
00454 const std::vector<G4VEnergyLossProcess*>& v =
00455 theManager->GetEnergyLossProcessVector();
00456 std::vector<G4VEnergyLossProcess*>::const_iterator itr;
00457 for(itr = v.begin(); itr != v.end(); ++itr) {
00458 G4VEnergyLossProcess* p = *itr;
00459 if(p) {
00460 if (p->GetProcessName() == name) {
00461 p->ActivateForcedInteraction(length,region,flag);
00462 }
00463 }
00464 }
00465 const std::vector<G4VEmProcess*>& w =
00466 theManager->GetEmProcessVector();
00467 std::vector<G4VEmProcess*>::const_iterator itp;
00468 for(itp = w.begin(); itp != w.end(); itp++) {
00469 G4VEmProcess* q = *itp;
00470 if(q) {
00471 if (q->GetProcessName() == name) {
00472 q->ActivateForcedInteraction(length,region,flag);
00473 }
00474 }
00475 }
00476 }
00477
00478
00479
00480 void
00481 G4EmProcessOptions::ActivateSecondaryBiasing(const G4String& name,
00482 const G4String& region,
00483 G4double factor,
00484 G4double energyLimit)
00485 {
00486 if(0.0 > factor) { return; }
00487 const std::vector<G4VEnergyLossProcess*>& v =
00488 theManager->GetEnergyLossProcessVector();
00489 std::vector<G4VEnergyLossProcess*>::const_iterator itr;
00490 for(itr = v.begin(); itr != v.end(); ++itr) {
00491 G4VEnergyLossProcess* p = *itr;
00492 if(p) {
00493 if (p->GetProcessName() == name) {
00494 p->ActivateSecondaryBiasing(region, factor, energyLimit);
00495 }
00496 }
00497 }
00498 }
00499
00500
00501
00502 void
00503 G4EmProcessOptions::ActivateSecondaryBiasingForGamma(const G4String& name,
00504 const G4String& region,
00505 G4double factor,
00506 G4double energyLimit)
00507 {
00508 if(0.0 > factor) { return; }
00509 const std::vector<G4VEmProcess*>& v =
00510 theManager->GetEmProcessVector();
00511 std::vector<G4VEmProcess*>::const_iterator itr;
00512 for(itr = v.begin(); itr != v.end(); ++itr) {
00513 G4VEmProcess* p = *itr;
00514 if(p) {
00515 if (p->GetProcessName() == name) {
00516 p->ActivateSecondaryBiasing(region, factor, energyLimit);
00517 }
00518 }
00519 }
00520 }
00521
00522
00523