Working with the G4Navigator class in user code

Hi,

I’m trying to set up an electric field based on a CAD file. The CAD file defines a set of cells. Each cell is associated to an E-field vector. When Geant4 asks for the value of the E-field in a given point, I need to figure out which cell contains the point and return the corresponding value. I’m hoping to achieve this by using a G4Navigator, but I want to make sure that my setup won’t cause any unexpected problems since the G4Navigator class isn’t normally supposed to be touched by the user.

I have a class MyElectricField which inherits from G4ElectricField. I could read the geometry and create a G4Navigator in the constructor.

void MyElectricField::MyElectricField() {
	fieldNavigator = new G4Navigator(); //fieldNavigator is a member variable of MyElectricField
	fieldNavigator->SetWorldVolume(worldPhysical); //worldPhysical has been populated with the CAD geometry
}

Then I could ask the G4Navigator to locate the desired volume each time GetFieldValue is called.

void MyElectricField::GetFieldValue(const G4double Point[4], G4double *field) const {
	G4VPhysicalVolume* volume = fieldNavigator->LocateGlobalPointAndSetup(G4ThreeVector(Point[0], Point[1], Point[2]));
	G4ThreeVector eField = eFieldVector.at(volume->GetCopyNo()); //The E-field in each cell is stored in a vector which is accessed by copyNo 
}
  1. Would this work as expected?

  2. Is there a better way to achieve what I want?

  3. Looking in the G4Navigator src, some functions (like LocateGlobalPointAndSetup) note that the geometry must be closed. How is this achieved? Do I need to involve the G4GeometryManager somehow to close and/or optimize the geometry?

Thanks in advance!

//
Anders

  1. This does work as expected. We’re using it in our simulation for a very similar purpose (to get from global points to volumes, and then to the stuff associated with that volume).

  2. G4Navigator is the right tool for the job. You’ve made it a data member of your electric field, which is good. The field object is already thread local.

  3. Unless you’re doing some peculiar diagnostics, GetFieldValue() will only be called during the event loop. The geometry is closed as part of /run/initialize, so what you have there is fine.

Hi Anders, (Mike)

I confirm Mike’s reply.

You can create a separate ‘ghost’ geometry that does not share volumes with the ‘mass’ geometry in which the setup’s volumes are defined, and then create one navigator in each thread with which to locate points. It’s volumes will be optimised for location/navigation, as all volumes registered in G4LogicalVolumeStore are scanned and optimised, not only the ones in the ‘mass’ geometry.

Please note that the GetFieldValue() could be called for positions which are outside the world volume - or the expected volume of the field. It will only check where they are located after the integration.

One caveat: in case you are sharing the ‘mass’ geometry with the field, then your geometry setup must not contain replicas, parameterised volumes or divisions.

Hi,

Thank you for your answers! Everything seems to work. The E-field geometry is indeed separate from the mass geometry. If LocateGlobalPointAndSetup returns zero, I handle it by simply returning a field of zero.

Here are some implementation details which may be useful to someone in the future. The DetectorConstruction keeps track of the geometry in which the E-field exists, and has a thread-local field manager.

class MyDetectorConstruction : public G4VUserDetectorConstruction {
	G4VPhysicalVolume* eFieldPhysical; //This is the mother physical volume of the eField geometry
	G4Cache<MyElectricFieldManager*> eFieldManager; //The ElectricField is a member variable of this class. I suppose making it a (G4Cache'd) member of MyDetectorConstruction would also have worked.
}

In ConstructSDandField, I read the CAD geometry and pass it to the ElectricField. I make sure to only build the geometry once. The mutex is there to keep many worker threads from building geometries at the same time, but I’m not sure it’s required as the master thread seems to call ConstructSDandField before any worker threads are even created.

void MyDetectorConstruction::ConstructSDandField() {
	MyElectricFieldManager* fieldManager = new MyElectricFieldManager();
	G4AutoLock lock(&MyDetectorConstructionMutex); //Not sure that this mutex is needed. It looks like the master thread calls ConstructSDandField before any worker threads do.
	if(!eFieldPhysical) eFieldPhysical = fieldManager->getField()->construct(); //We only need to build the CAD geometry once
	lock.unlock();
	
	fieldManager->getField()->setGeometry(eFieldPhysical); //The field will create a G4Navigator to navigate this geometry
	G4AutoDelete::Register(fieldManager); //The kernel will delete this object
	eFieldManager.Put(fieldManager);
}

For the sake of my curiosity, here are two more questions.

  1. Is the mutex needed, or can I always rely on the master thread to set eFieldPhysical before any worker threads do?

  2. I looked at example basic/B5, and they declared the eFieldManager as static G4ThreadLocal rather than using G4Cache. Does it matter in this case? My understanding is that G4ThreadLocal is for variables that should be shared among many instances of MyDetectorConstruction whereas G4Cache is for variables that should not be shared. I suppose only one instance of MyDetectorConstruction exists, so perhaps it is a moot point.

//
Anders