Message: Re: Uniform magnetic field using the G4GlobalMagFieldMessenger Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Uniform magnetic field using the G4GlobalMagFieldMessenger 

Forum: Fields: Magnetic and Otherwise
Re: None Uniform magnetic field using the G4GlobalMagFieldMessenger
Date: 04 Apr, 2018
From: Eric Simiele <Eric Simiele>

Hello,

I would recommend having a look at extended example field01. Essentially, you should write your own magnetic field class similar to F01FieldSetup where you can specify all of the parameters related to the magnetic field (e.g., chord finder, min step in field, max step in field, miss distance, etc.). You would instantiate this class from the ConstructSDandField method of your DetectorConstruction class. Furthermore, you can write a messenger class for this magnetic field class so you can change the various field-specific parameters at run time.

Good luck!

Eric

On Mon, 19 Feb 2018 19:06:27 GMT, heid wrote:

> Hi
> 
> I am trying to simulate (with Geant4 10.4 and multithreading) a uniform
> magnetic field using the G4GlobalMagFieldMessenger in my
> ConstructSDandField method of DetectorConstruction.cc (see code below)
> and the UIcommand /globalField/setValue 0 1.5 0 tesla in my input file.
> My problem is that I would like to modify the accuracy of the fMinStep
> transport parameter, but I am having difficulty with this.
> 
> Looking at the Geant4 TestEM examples, I cannot work out how to access
> and manipulate the G4ChordFinder set in the G4GlobalMagFieldMessenger
> within my ConstructSDandField method of DetectorConstruction.cc in order
> to change the value of fMinStep.
> 
> I am new to Geant4 and am hoping that someone could please steer me in
> the right direction.
> 
> --------------------------------------------------------------------------------
> G4ThreadLocal G4GlobalMagFieldMessenger* DetectorConstruction::fMagFieldMsg = 0;
> G4ThreadLocal G4UniformMagField* DetectorConstruction::MagField = 0;        
> G4ThreadLocal G4FieldManager* DetectorConstruction::fieldMgr = 0;
> void DetectorConstruction::ConstructSDandField()
> {
>   // Only called by Master thread
>   std::ios::fmtflags mode = G4cout.flags();
>   G4cout.setf(std::ios::fixed,std::ios::floatfield);
> 
>   //Define a global uniform magnetic field
>   if (!fMagFieldMsg)
>     { 
>       //--- Set value of accuracy parameters ------------------------------
>       G4double minEps        = 5.0e-5;    // Minimum & value for smallest steps
>       G4double maxEps        = 5.0e-5;    // Maximum & value for largest steps
>       G4double delta_onestep = 1.0*micrometer; // 1.0 micrometer
>       G4double miss_distance = 1.0*micrometer;
>       G4double delta_intsect = 1.0*micrometer;
>       G4double fMinStep      = 1.0*nm;    // minimal step of 1nm
>       G4double fMaxStep      = 1.0e12;   
>       G4MagIntegratorStepper* stepper = 0;//G4 default stepper
> 
>       //--- Create global magnetic field messenger. 
>       //Uniform magnetic field AUTO created when field is non-zero.
>       G4ThreeVector fieldValue = G4ThreeVector();
>       fMagFieldMsg = new G4GlobalMagFieldMessenger(fieldValue);
>       fMagFieldMsg->SetVerboseLevel(1);
> 
>       //Get value of Magnetic field from global magnetic field messenger.
>       G4ThreeVector magField =  fMagFieldMsg->GetFieldValue();
> 
>       G4FieldManager *fieldMgr;
>       G4TransportationManager *transportMgr =
>                       G4TransportationManager::GetTransportationManager();
>       fieldMgr = transportMgr->GetFieldManager();
>       //fieldMgr->SetDetectorField(magField);
>       //fieldMgr->CreateChordFinder(magField)
> 
>       transportMgr->GetPropagatorInField()->SetLargestAcceptableStep(fMaxStep);
> 
>       fieldMgr->SetMinimumEpsilonStep( minEps );
>       fieldMgr->SetMaximumEpsilonStep( maxEps );
>       fieldMgr->SetDeltaOneStep( delta_onestep );  // 0.5 micrometer
>       fieldMgr->GetChordFinder()->SetDeltaChord( miss_distance );
>       fieldMgr->SetDeltaIntersection( delta_intsect );
>       //fieldMgr->SetChordFinder(fChordFinder);
> 
>       //Register the field messenger for deleting
>       G4AutoDelete::Register(fMagFieldMsg);
>    }
> }
> --------------------------------------------------------------------------------
> 

Inline Depth:
 1 1
 All All
Outline Depth:
 1 1
 2 2
 All All
Add message: (add)

1 None: Re: Uniform magnetic field using the G4GlobalMagFieldMessenger   (heid - 05 Apr, 2018)
 Add Message Add Message
to: "Re: Uniform magnetic field using the G4GlobalMagFieldMessenger"

 Subscribe Subscribe

This site runs SLAC HyperNews version 1.11-slac-98, derived from the original HyperNews


[ Geant 4 Home | Geant 4 HyperNews | Search | Request New Forum | Feedback ]