#include <G4Decay.hh>
Inheritance diagram for G4Decay:
Definition at line 63 of file G4Decay.hh.
G4Decay::G4Decay | ( | const G4String & | processName = "Decay" |
) |
Definition at line 63 of file G4Decay.cc.
References DECAY, fParticleChangeForDecay, G4cout, G4endl, GetVerboseLevel(), G4VProcess::pParticleChange, and G4VProcess::SetProcessSubType().
00064 :G4VRestDiscreteProcess(processName, fDecay), 00065 verboseLevel(1), 00066 HighestValue(20.0), 00067 fRemainderLifeTime(-1.0), 00068 pExtDecayer(0) 00069 { 00070 // set Process Sub Type 00071 SetProcessSubType(static_cast<int>(DECAY)); 00072 00073 #ifdef G4VERBOSE 00074 if (GetVerboseLevel()>1) { 00075 G4cout << "G4Decay constructor " << " Name:" << processName << G4endl; 00076 } 00077 #endif 00078 00079 pParticleChange = &fParticleChangeForDecay; 00080 }
G4Decay::~G4Decay | ( | ) | [virtual] |
Definition at line 82 of file G4Decay.cc.
References pExtDecayer.
00083 { 00084 if (pExtDecayer) { 00085 delete pExtDecayer; 00086 } 00087 }
G4VParticleChange * G4Decay::AtRestDoIt | ( | const G4Track & | aTrack, | |
const G4Step & | aStep | |||
) | [inline, virtual] |
Reimplemented from G4VRestDiscreteProcess.
Definition at line 191 of file G4Decay.hh.
References DecayIt().
00195 { 00196 return DecayIt(aTrack, aStep); 00197 }
G4double G4Decay::AtRestGetPhysicalInteractionLength | ( | const G4Track & | track, | |
G4ForceCondition * | condition | |||
) | [virtual] |
Reimplemented from G4VRestDiscreteProcess.
Definition at line 438 of file G4Decay.cc.
References DBL_MIN, fRemainderLifeTime, G4Track::GetDynamicParticle(), GetMeanLifeTime(), G4DynamicParticle::GetPreAssignedDecayProperTime(), G4Track::GetProperTime(), NotForced, and G4VProcess::theNumberOfInteractionLengthLeft.
00442 { 00443 // condition is set to "Not Forced" 00444 *condition = NotForced; 00445 00446 G4double pTime = track.GetDynamicParticle()->GetPreAssignedDecayProperTime(); 00447 if (pTime >= 0.) { 00448 fRemainderLifeTime = pTime - track.GetProperTime(); 00449 if (fRemainderLifeTime <= 0.0) fRemainderLifeTime = DBL_MIN; 00450 } else { 00451 fRemainderLifeTime = 00452 theNumberOfInteractionLengthLeft * GetMeanLifeTime(track, condition); 00453 } 00454 return fRemainderLifeTime; 00455 }
void G4Decay::BuildPhysicsTable | ( | const G4ParticleDefinition & | ) | [virtual] |
void G4Decay::DaughterPolarization | ( | const G4Track & | aTrack, | |
G4DecayProducts * | products | |||
) | [protected, virtual] |
Reimplemented in G4PionDecayMakeSpin.
Definition at line 348 of file G4Decay.cc.
Referenced by DecayIt().
G4VParticleChange * G4Decay::DecayIt | ( | const G4Track & | aTrack, | |
const G4Step & | aStep | |||
) | [protected, virtual] |
Reimplemented in G4DecayWithSpin.
Definition at line 180 of file G4Decay.cc.
References G4VParticleChange::AddSecondary(), G4DecayProducts::Boost(), G4VProcess::ClearNumberOfInteractionLengthLeft(), DaughterPolarization(), G4VDecayChannel::DecayIt(), G4DecayProducts::DumpInfo(), G4DecayProducts::entries(), FatalException, fParticleChangeForDecay, fRemainderLifeTime, fStopAndKill, fStopButAlive, G4cout, G4endl, G4Exception(), G4ParticleDefinition::GetDecayTable(), G4DynamicParticle::GetDefinition(), G4Track::GetDynamicParticle(), G4Track::GetGlobalTime(), G4DynamicParticle::GetKineticEnergy(), G4Track::GetLocalTime(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetPDGStable(), G4Track::GetPosition(), G4DynamicParticle::GetPreAssignedDecayProducts(), G4DynamicParticle::GetTotalEnergy(), G4Track::GetTouchableHandle(), G4Track::GetTrackStatus(), G4VDecayChannel::GetVerboseLevel(), GetVerboseLevel(), G4VExtDecayer::ImportDecayProducts(), G4ParticleChangeForDecay::Initialize(), G4DecayProducts::IsChecked(), JustWarning, ns, pExtDecayer, G4DecayProducts::PopProducts(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForDecay::ProposeLocalTime(), G4VParticleChange::ProposeTrackStatus(), G4DecayTable::SelectADecayChannel(), G4Track::SetGoodForTrackingFlag(), G4VParticleChange::SetNumberOfSecondaries(), G4Track::SetTouchableHandle(), and G4VDecayChannel::SetVerboseLevel().
Referenced by AtRestDoIt(), G4DecayWithSpin::DecayIt(), and PostStepDoIt().
00181 { 00182 // The DecayIt() method returns by pointer a particle-change object. 00183 // Units are expressed in GEANT4 internal units. 00184 00185 // Initialize ParticleChange 00186 // all members of G4VParticleChange are set to equal to 00187 // corresponding member in G4Track 00188 fParticleChangeForDecay.Initialize(aTrack); 00189 00190 // get particle 00191 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle(); 00192 const G4ParticleDefinition* aParticleDef = aParticle->GetDefinition(); 00193 00194 // check if the particle is stable 00195 if (aParticleDef->GetPDGStable()) return &fParticleChangeForDecay ; 00196 00197 00198 //check if thePreAssignedDecayProducts exists 00199 const G4DecayProducts* o_products = (aParticle->GetPreAssignedDecayProducts()); 00200 G4bool isPreAssigned = (o_products != 0); 00201 G4DecayProducts* products = 0; 00202 00203 // decay table 00204 G4DecayTable *decaytable = aParticleDef->GetDecayTable(); 00205 00206 // check if external decayer exists 00207 G4bool isExtDecayer = (decaytable == 0) && (pExtDecayer !=0); 00208 00209 // Error due to NO Decay Table 00210 if ( (decaytable == 0) && !isExtDecayer &&!isPreAssigned ){ 00211 if (GetVerboseLevel()>0) { 00212 G4cout << "G4Decay::DoIt : decay table not defined for "; 00213 G4cout << aParticle->GetDefinition()->GetParticleName()<< G4endl; 00214 } 00215 G4Exception( "G4Decay::DecayIt ", 00216 "DECAY101",JustWarning, 00217 "Decay table is not defined"); 00218 00219 fParticleChangeForDecay.SetNumberOfSecondaries(0); 00220 // Kill the parent particle 00221 fParticleChangeForDecay.ProposeTrackStatus( fStopAndKill ) ; 00222 fParticleChangeForDecay.ProposeLocalEnergyDeposit(0.0); 00223 00224 ClearNumberOfInteractionLengthLeft(); 00225 return &fParticleChangeForDecay ; 00226 } 00227 00228 if (isPreAssigned) { 00229 // copy decay products 00230 products = new G4DecayProducts(*o_products); 00231 } else if ( isExtDecayer ) { 00232 // decay according to external decayer 00233 products = pExtDecayer->ImportDecayProducts(aTrack); 00234 } else { 00235 // decay acoording to decay table 00236 // choose a decay channel 00237 G4VDecayChannel *decaychannel = decaytable->SelectADecayChannel(); 00238 if (decaychannel == 0 ){ 00239 // decay channel not found 00240 G4Exception("G4Decay::DoIt", "DECAY003", FatalException, 00241 " can not determine decay channel "); 00242 } else { 00243 // execute DecayIt() 00244 #ifdef G4VERBOSE 00245 G4int temp = decaychannel->GetVerboseLevel(); 00246 if (GetVerboseLevel()>1) { 00247 G4cout << "G4Decay::DoIt : selected decay channel addr:" << decaychannel <<G4endl; 00248 decaychannel->SetVerboseLevel(GetVerboseLevel()); 00249 } 00250 #endif 00251 products = decaychannel->DecayIt(aParticle->GetMass()); 00252 #ifdef G4VERBOSE 00253 if (GetVerboseLevel()>1) { 00254 decaychannel->SetVerboseLevel(temp); 00255 } 00256 #endif 00257 #ifdef G4VERBOSE 00258 if (GetVerboseLevel()>2) { 00259 if (! products->IsChecked() ) products->DumpInfo(); 00260 } 00261 #endif 00262 } 00263 } 00264 00265 // get parent particle information ................................... 00266 G4double ParentEnergy = aParticle->GetTotalEnergy(); 00267 G4double ParentMass = aParticle->GetMass(); 00268 if (ParentEnergy < ParentMass) { 00269 if (GetVerboseLevel()>0) { 00270 G4cout << "G4Decay::DoIt : Total Energy is less than its mass" << G4endl; 00271 G4cout << " Particle: " << aParticle->GetDefinition()->GetParticleName(); 00272 G4cout << " Energy:" << ParentEnergy/MeV << "[MeV]"; 00273 G4cout << " Mass:" << ParentMass/MeV << "[MeV]"; 00274 G4cout << G4endl; 00275 } 00276 G4Exception( "G4Decay::DecayIt ", 00277 "DECAY102",JustWarning, 00278 "Total Energy is less than its mass"); 00279 ParentEnergy = ParentMass; 00280 } 00281 00282 G4ThreeVector ParentDirection(aParticle->GetMomentumDirection()); 00283 00284 //boost all decay products to laboratory frame 00285 G4double energyDeposit = 0.0; 00286 G4double finalGlobalTime = aTrack.GetGlobalTime(); 00287 G4double finalLocalTime = aTrack.GetLocalTime(); 00288 if (aTrack.GetTrackStatus() == fStopButAlive ){ 00289 // AtRest case 00290 finalGlobalTime += fRemainderLifeTime; 00291 finalLocalTime += fRemainderLifeTime; 00292 energyDeposit += aParticle->GetKineticEnergy(); 00293 if (isPreAssigned) products->Boost( ParentEnergy, ParentDirection); 00294 } else { 00295 // PostStep case 00296 if (!isExtDecayer) products->Boost( ParentEnergy, ParentDirection); 00297 } 00298 00299 // set polarization for daughter particles 00300 DaughterPolarization(aTrack, products); 00301 00302 00303 //add products in fParticleChangeForDecay 00304 G4int numberOfSecondaries = products->entries(); 00305 fParticleChangeForDecay.SetNumberOfSecondaries(numberOfSecondaries); 00306 #ifdef G4VERBOSE 00307 if (GetVerboseLevel()>1) { 00308 G4cout << "G4Decay::DoIt : Decay vertex :"; 00309 G4cout << " Time: " << finalGlobalTime/ns << "[ns]"; 00310 G4cout << " X:" << (aTrack.GetPosition()).x() /cm << "[cm]"; 00311 G4cout << " Y:" << (aTrack.GetPosition()).y() /cm << "[cm]"; 00312 G4cout << " Z:" << (aTrack.GetPosition()).z() /cm << "[cm]"; 00313 G4cout << G4endl; 00314 G4cout << "G4Decay::DoIt : decay products in Lab. Frame" << G4endl; 00315 products->DumpInfo(); 00316 } 00317 #endif 00318 G4int index; 00319 G4ThreeVector currentPosition; 00320 const G4TouchableHandle thand = aTrack.GetTouchableHandle(); 00321 for (index=0; index < numberOfSecondaries; index++) 00322 { 00323 // get current position of the track 00324 currentPosition = aTrack.GetPosition(); 00325 // create a new track object 00326 G4Track* secondary = new G4Track( products->PopProducts(), 00327 finalGlobalTime , 00328 currentPosition ); 00329 // switch on good for tracking flag 00330 secondary->SetGoodForTrackingFlag(); 00331 secondary->SetTouchableHandle(thand); 00332 // add the secondary track in the List 00333 fParticleChangeForDecay.AddSecondary(secondary); 00334 } 00335 delete products; 00336 00337 // Kill the parent particle 00338 fParticleChangeForDecay.ProposeTrackStatus( fStopAndKill ) ; 00339 fParticleChangeForDecay.ProposeLocalEnergyDeposit(energyDeposit); 00340 fParticleChangeForDecay.ProposeLocalTime( finalLocalTime ); 00341 00342 // Clear NumberOfInteractionLengthLeft 00343 ClearNumberOfInteractionLengthLeft(); 00344 00345 return &fParticleChangeForDecay ; 00346 }
void G4Decay::EndTracking | ( | ) | [virtual] |
Reimplemented from G4VProcess.
Definition at line 362 of file G4Decay.cc.
References G4VProcess::ClearNumberOfInteractionLengthLeft(), and G4VProcess::currentInteractionLength.
00363 { 00364 // Clear NumberOfInteractionLengthLeft 00365 ClearNumberOfInteractionLengthLeft(); 00366 00367 currentInteractionLength = -1.0; 00368 }
const G4VExtDecayer * G4Decay::GetExtDecayer | ( | ) | const [inline] |
Definition at line 209 of file G4Decay.hh.
References pExtDecayer.
00210 { 00211 return pExtDecayer; 00212 }
G4double G4Decay::GetMeanFreePath | ( | const G4Track & | aTrack, | |
G4double | previousStepSize, | |||
G4ForceCondition * | condition | |||
) | [protected, virtual] |
Implements G4VRestDiscreteProcess.
Definition at line 129 of file G4Decay.cc.
References DBL_MAX, DBL_MIN, G4cout, G4endl, G4DynamicParticle::GetDefinition(), G4Track::GetDynamicParticle(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMass(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetPDGLifeTime(), G4ParticleDefinition::GetPDGStable(), G4DynamicParticle::GetTotalMomentum(), GetVerboseLevel(), and HighestValue.
Referenced by PostStepGetPhysicalInteractionLength().
00130 { 00131 // get particle 00132 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle(); 00133 const G4ParticleDefinition* aParticleDef = aParticle->GetDefinition(); 00134 G4double aMass = aParticle->GetMass(); 00135 G4double aLife = aParticleDef->GetPDGLifeTime(); 00136 00137 00138 // returns the mean free path in GEANT4 internal units 00139 G4double pathlength; 00140 G4double aCtau = c_light * aLife; 00141 00142 // check if the particle is stable? 00143 if (aParticleDef->GetPDGStable()) { 00144 pathlength = DBL_MAX; 00145 00146 //check if the particle has very short life time ? 00147 } else if (aCtau < DBL_MIN) { 00148 pathlength = DBL_MIN; 00149 00150 } else { 00151 //calculate the mean free path 00152 // by using normalized kinetic energy (= Ekin/mass) 00153 G4double rKineticEnergy = aParticle->GetKineticEnergy()/aMass; 00154 if ( rKineticEnergy > HighestValue) { 00155 // gamma >> 1 00156 pathlength = ( rKineticEnergy + 1.0)* aCtau; 00157 } else if ( rKineticEnergy < DBL_MIN ) { 00158 // too slow particle 00159 #ifdef G4VERBOSE 00160 if (GetVerboseLevel()>1) { 00161 G4cout << "G4Decay::GetMeanFreePath() !!particle stops!!"; 00162 G4cout << aParticleDef->GetParticleName() << G4endl; 00163 G4cout << "KineticEnergy:" << aParticle->GetKineticEnergy()/GeV <<"[GeV]"; 00164 } 00165 #endif 00166 pathlength = DBL_MIN; 00167 } else { 00168 // beta <1 00169 pathlength = (aParticle->GetTotalMomentum())/aMass*aCtau ; 00170 } 00171 } 00172 return pathlength; 00173 }
G4double G4Decay::GetMeanLifeTime | ( | const G4Track & | aTrack, | |
G4ForceCondition * | condition | |||
) | [protected, virtual] |
Implements G4VRestDiscreteProcess.
Definition at line 101 of file G4Decay.cc.
References DBL_MAX, G4cout, G4endl, G4DynamicParticle::GetDefinition(), G4Track::GetDynamicParticle(), G4ParticleDefinition::GetPDGLifeTime(), G4ParticleDefinition::GetPDGStable(), GetVerboseLevel(), and ns.
Referenced by AtRestGetPhysicalInteractionLength().
00103 { 00104 // returns the mean free path in GEANT4 internal units 00105 G4double meanlife; 00106 00107 // get particle 00108 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle(); 00109 const G4ParticleDefinition* aParticleDef = aParticle->GetDefinition(); 00110 G4double aLife = aParticleDef->GetPDGLifeTime(); 00111 00112 // check if the particle is stable? 00113 if (aParticleDef->GetPDGStable()) { 00114 meanlife = DBL_MAX; 00115 00116 } else { 00117 meanlife = aLife; 00118 } 00119 00120 #ifdef G4VERBOSE 00121 if (GetVerboseLevel()>1) { 00122 G4cout << "mean life time: "<< meanlife/ns << "[ns]" << G4endl; 00123 } 00124 #endif 00125 00126 return meanlife; 00127 }
G4double G4Decay::GetRemainderLifeTime | ( | ) | const [inline] |
Definition at line 215 of file G4Decay.hh.
References fRemainderLifeTime.
00216 { 00217 return fRemainderLifeTime; 00218 }
G4int G4Decay::GetVerboseLevel | ( | ) | const [inline] |
Reimplemented from G4VProcess.
Definition at line 188 of file G4Decay.hh.
References verboseLevel.
Referenced by DecayIt(), G4Decay(), GetMeanFreePath(), and GetMeanLifeTime().
00188 { return verboseLevel; }
G4bool G4Decay::IsApplicable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VProcess.
Definition at line 89 of file G4Decay.cc.
References G4ParticleDefinition::GetPDGLifeTime(), and G4ParticleDefinition::GetPDGMass().
Referenced by TLBE< T >::ConstructGeneral(), and G4DecayPhysics::ConstructProcess().
00090 { 00091 // check if the particle is stable? 00092 if (aParticleType.GetPDGLifeTime() <0.0) { 00093 return false; 00094 } else if (aParticleType.GetPDGMass() <= 0.0*MeV) { 00095 return false; 00096 } else { 00097 return true; 00098 } 00099 }
G4VParticleChange * G4Decay::PostStepDoIt | ( | const G4Track & | aTrack, | |
const G4Step & | aStep | |||
) | [inline, virtual] |
Reimplemented from G4VRestDiscreteProcess.
Definition at line 200 of file G4Decay.hh.
References DecayIt().
00204 { 00205 return DecayIt(aTrack, aStep); 00206 }
G4double G4Decay::PostStepGetPhysicalInteractionLength | ( | const G4Track & | track, | |
G4double | previousStepSize, | |||
G4ForceCondition * | condition | |||
) | [virtual] |
Reimplemented from G4VRestDiscreteProcess.
Definition at line 371 of file G4Decay.cc.
References G4VProcess::currentInteractionLength, DBL_MAX, DBL_MIN, G4DynamicParticle::DumpInfo(), fRemainderLifeTime, G4cout, G4endl, G4DynamicParticle::GetDefinition(), G4Track::GetDynamicParticle(), G4DynamicParticle::GetMass(), G4Track::GetMaterial(), GetMeanFreePath(), G4Material::GetName(), G4ParticleDefinition::GetPDGLifeTime(), G4DynamicParticle::GetPreAssignedDecayProperTime(), G4Track::GetProperTime(), G4DynamicParticle::GetTotalMomentum(), NotForced, G4VProcess::SubtractNumberOfInteractionLengthLeft(), G4VProcess::theNumberOfInteractionLengthLeft, and verboseLevel.
00376 { 00377 00378 // condition is set to "Not Forced" 00379 *condition = NotForced; 00380 00381 // pre-assigned Decay time 00382 G4double pTime = track.GetDynamicParticle()->GetPreAssignedDecayProperTime(); 00383 G4double aLife = track.GetDynamicParticle()->GetDefinition()->GetPDGLifeTime(); 00384 00385 if (pTime < 0.) { 00386 // normal case 00387 if ( previousStepSize > 0.0){ 00388 // subtract NumberOfInteractionLengthLeft 00389 SubtractNumberOfInteractionLengthLeft(previousStepSize); 00390 if(theNumberOfInteractionLengthLeft<0.){ 00391 theNumberOfInteractionLengthLeft=perMillion; 00392 } 00393 fRemainderLifeTime = theNumberOfInteractionLengthLeft*aLife; 00394 } 00395 // get mean free path 00396 currentInteractionLength = GetMeanFreePath(track, previousStepSize, condition); 00397 00398 #ifdef G4VERBOSE 00399 if ((currentInteractionLength <=0.0) || (verboseLevel>2)){ 00400 G4cout << "G4Decay::PostStepGetPhysicalInteractionLength " << G4endl; 00401 track.GetDynamicParticle()->DumpInfo(); 00402 G4cout << " in Material " << track.GetMaterial()->GetName() <<G4endl; 00403 G4cout << "MeanFreePath = " << currentInteractionLength/cm << "[cm]" <<G4endl; 00404 } 00405 #endif 00406 00407 G4double value; 00408 if (currentInteractionLength <DBL_MAX) { 00409 value = theNumberOfInteractionLengthLeft * currentInteractionLength; 00410 } else { 00411 value = DBL_MAX; 00412 } 00413 00414 return value; 00415 00416 } else { 00417 //pre-assigned Decay time case 00418 // reminder proper time 00419 fRemainderLifeTime = pTime - track.GetProperTime(); 00420 if (fRemainderLifeTime <= 0.0) fRemainderLifeTime = DBL_MIN; 00421 00422 G4double rvalue=0.0; 00423 // use pre-assigned Decay time to determine PIL 00424 if (aLife>0.0) { 00425 // ordinary particle 00426 rvalue = (fRemainderLifeTime/aLife)*GetMeanFreePath(track, previousStepSize, condition); 00427 } else { 00428 // shortlived particle 00429 rvalue = c_light * fRemainderLifeTime; 00430 // by using normalized kinetic energy (= Ekin/mass) 00431 G4double aMass = track.GetDynamicParticle()->GetMass(); 00432 rvalue *= track.GetDynamicParticle()->GetTotalMomentum()/aMass; 00433 } 00434 return rvalue; 00435 } 00436 }
void G4Decay::SetExtDecayer | ( | G4VExtDecayer * | ) |
Definition at line 458 of file G4Decay.cc.
References DECAY_External, pExtDecayer, and G4VProcess::SetProcessSubType().
00459 { 00460 pExtDecayer = val; 00461 00462 // set Process Sub Type 00463 if ( pExtDecayer !=0 ) { 00464 SetProcessSubType(static_cast<int>(DECAY_External)); 00465 } 00466 }
void G4Decay::SetVerboseLevel | ( | G4int | value | ) | [inline] |
Reimplemented from G4VProcess.
Definition at line 185 of file G4Decay.hh.
References verboseLevel.
00185 { verboseLevel = value; }
void G4Decay::StartTracking | ( | G4Track * | ) | [virtual] |
Reimplemented from G4VProcess.
Definition at line 354 of file G4Decay.cc.
References G4VProcess::currentInteractionLength, fRemainderLifeTime, and G4VProcess::ResetNumberOfInteractionLengthLeft().
00355 { 00356 currentInteractionLength = -1.0; 00357 ResetNumberOfInteractionLengthLeft(); 00358 00359 fRemainderLifeTime = -1.0; 00360 }
G4double G4Decay::fRemainderLifeTime [protected] |
Definition at line 175 of file G4Decay.hh.
Referenced by AtRestGetPhysicalInteractionLength(), G4DecayWithSpin::DecayIt(), DecayIt(), GetRemainderLifeTime(), PostStepGetPhysicalInteractionLength(), and StartTracking().
const G4double G4Decay::HighestValue [protected] |
G4VExtDecayer* G4Decay::pExtDecayer [protected] |
Definition at line 181 of file G4Decay.hh.
Referenced by DecayIt(), GetExtDecayer(), SetExtDecayer(), and ~G4Decay().
G4int G4Decay::verboseLevel [protected] |
Reimplemented from G4VProcess.
Definition at line 164 of file G4Decay.hh.
Referenced by GetVerboseLevel(), PostStepGetPhysicalInteractionLength(), and SetVerboseLevel().