Generate particles from a volume

Hi,

very general question. I have logical and physical volumes defining the geometry and I want to generate events from some of those volumes without re-defining particles generators. Is that possible?

The problem of course must be common: having materials with some known radiation activity and want to simulate that knowing an energy spectrum.

Can I associate to a specific volume some particle generator?

Thanks,
Paolo

Unfortunately, this is not available directly in the toolkit. G4GeneralParticleSource (GPS) supports defining shapes, but not specifying a named volume in your geometry. There might be instances in the extended or advanced examples that show how to do what you want.

We have implemented this in our experiment’s simulation framework for exactly the reason you specify – measuring the effect of radiocontaminants on our detectors. To do it, you need to write your own generator (which you could make either a subclass of, or owner of an instance of, G4ParticleGun).

  • Make the volume name configurable, and you can look it up in the appropriate G4*VolumeStore.
  • Use G4VSolid::BoundingLimits() to get a simple region enclosing the volume, and use an accept-reject loop with G4VSolid::Inside() to select (local coordinate system) points in the bulk of the volume.
  • Be sure you transform from local to global coordinates for the generator.

Thanks for your help.
So If I understand correctly, I assume you loop within the boundaries coordinate of a selected volume (with some stepping) to generate a point-like particle gun in random directions?

Thanks!

Not exactly, as what you describe won’t give you a uniform distribution. Here’s some pseudo-code for generating uniform points in a volume.

Let (xlo,ylo,zlo) and (xhi,yhi,zhi) define a box that fully encloses your volume (which might be a sphere, or have holes in it, or whatever). Let pvol be a pointer to the specific placement volume you want to “contaminate.”

G4VSolid solid = pvol->GetLogicalVolume()->GetSolid();
G4ThreeVector point;
G4int maxtries=10000, itry=1;
do {
  point.set(xlo + G4UniformRand()*(xhi-xlo),
            ylo + G4UniformRand()*(yhi-ylo),
            zlo + G4UniformRand()*(zhi-zlo));
} while (!solid->Inside(point) && ++itry < maxtries);

if (itry == maxtries)
  G4cerr << "Unable to find a point inside your volume!" << G4endl;

If the loop doesn’t fail, then point will have a uniform random distribution within your volume.

Ok thanks for the clarification. So to confirm, for each point you instruct a G4ParticleGun with some random direction.

Yes; that’s how I would do it, if you’re generating gammas, or electrons, or whatever.

If you’re placing a radioisotope (like Si-32, or tritium, or something in a chip), then you can leave direction as (0,0,0), since an atom placement would be at rest.

Thanks, I will write a function to input any volume I want to use as a generator.

To be most general as possible I was trying to get the boundaries of the enclosure with

G4VisExtent *extent = solid->GetExtent();

but apparently I get something like this

invalid use of incomplete type ‘class G4VisExtent’ which is not clear to me how should be wrong.

Thanks!

“Incomplete type” is compiler shorthand for, “You forgot the .hh file, you git!” :slight_smile:

G4VisExtent might be useful, but it seems a bit heavyweight. For this purpose, you’ll be generating points in the local frame (and transforming them to world coordinates afterward). I think you could probably get away with

G4ThreeVector lo, hi;
solid->BoundingLimits(lo, hi);

Then use lo.x(), hi.y(), etc.

Thanks Michael, BoundingLimits did the job!

I am trying to access from my primary generator class the list of volumes defined in the detector construction class (and initialized by the runManager in the main). How can I do that without making the volumes global variables in some rather nasty way?

Of course the purpose of this, is to get any solid property I need before instructing the particle gun. Is there any kernel volume store?

Thanks again!

Yes, the geometry is kept in static (thread local) volume stores: G4SolidStore, G4LogicalVolumeStore, G4PhysicalVolumeStore, G4AssemblyStore, G4RegionStore. I’d recommend going in via the PV store to get the exactly volume (name and copy number) that you want, and then get to the Solid that belongs to that PV.

Many years ago, Pedro Arce (@arce) gave me a set of utilities from the GAMOS experiment that could take care of all this searching behind the scenes. We added the ability to get a full G4Touchable object for a volume, which then gives you the full global-to-local coordinate transforms you’ll need.

Let me grab our CDMS version of this tool and attach it here.

Here are the two classes we wrote, adapting the GAMOS tools.

DemoGeometryTools.cc (44.8 KB)
DemoGeometryTools.hh (6.8 KB)
DemoTouchable.cc (11.7 KB)
DemoTouchable.hh (4.4 KB)

The references to “DemoAbortJob” you’ll see are just a wrapper around calling G4Exception with FatalException. You can replace those with whatever error handling you want.

Thanks Mike, is all very much useful, will take great advantage out of them!

Speaking of radioactive sources, is there a way to do something like this
G4IonTable::GetIonTable()->GetIon("Co60");
instead of providing the pair (Z,A)?

thanks

Unfortunately not. I could imaging writing your own lookup table, but it won’t be very performant (string comparisons generally aren’t), parsing the G4-specified ion name string (“AbNN[e.eee]”) to get a simple Z,A,E triplet, but the name structure in general is sufficiently complicated that this probably isn’t worth it.

Considering the funny names in the NIST compounds list, I just assumed something could have been done for the isotopes, which have a less funky nomenclature.

If I provide a limited and specified list of names of isotopes (I do not have more of a dozen I am interested in) in the format AbNN will save any other user of my code (who is not me) the hassle of setting up Z,A in the configuration file.

Ah, but the NIST compounds are for the materials table, not for the particles list. You can’t do a “molecular beam” in Geant4, for example.

Yup! That’s sort of what I was imagining, and we do something similar in our simulation framework. If you’re writing a Messenger class for your generator anyway, providing support for selected ion names is straightforward.

yes sure, I was thinking just in terms of having nomenclatures like G4_CO60 for the G4IonTable.