BLCMDparticlefilter Class Reference

Inheritance diagram for BLCMDparticlefilter:

BLElement BLCommand

List of all members.


Detailed Description

class BLCMDparticlefilter - implements a ParticleFilter Each placement of this element generates a physical volume that will force the decay of any particle named in the "decay" argument, and will kill any particle named in its "kill" argument, or not named in its "keep" argument if that is non-empty.

Note that the standard Decay process is disabled for all particles named in the "decay" agument, so they will only decay when they enter a ParticleFilter volume.

Note also the G4VProcess of the ParticleFilterPlacement is only applied to particles listed in the decay argument.

There is a bug: the decay happens after the first step in the volume; so maxStep is arranged to ensure two steps occur in the volume. This appears to be in the way Geant4 handles physical volumes while tracking....

Public Member Functions

 BLCMDparticlefilter ()
 Default constructor.
virtual ~BLCMDparticlefilter ()
 Destructor.
BLElementclone ()
 clone()
 BLCMDparticlefilter (const BLCMDparticlefilter &r)
 Copy constructor.
G4String commandName ()
 commandName() returns "particlefilter".
int command (BLArgumentVector &argv, BLArgumentMap &namedArgs)
 command() implements the particlefilter command.
void defineNamedArgs ()
 defineNamedArgs() defines the named arguments for the command.
virtual void construct (G4RotationMatrix *relativeRotation, G4ThreeVector relativePosition, G4LogicalVolume *parent, G4String parentName, G4RotationMatrix *parentRotation, G4ThreeVector parentPosition)
 construct() will construct the particlefilter.
virtual G4double getLength ()
 getLength() returns this element's Length along the Z axis.
virtual G4double getWidth ()
 getWidth() returns this element's Width along the X axis.
virtual G4double getHeight ()
 getHeight() returns this element's height along the Y axis.
virtual G4bool isOK ()
 isOK() returns true.
void generatePoints (int npoints, std::vector< G4ThreeVector > &v)
 generatePoints() from BLElement
G4bool isOutside (G4ThreeVector &local, G4double tolerance)
 isOutside() from BLElement

Private Member Functions

void setup ()

Private Attributes

G4double radius
G4double innerRadius
G4double height
G4double width
G4double length
G4String material
G4String color
G4VSolid * solid
G4double maxStep
G4String decay
G4String kill
G4String keep
G4int nWait
G4int referenceWait
G4String require
G4int steppingVerbose
std::set< G4String > decaySet
std::set< G4String > killSet
std::set< G4String > keepSet

Static Private Attributes

static int instance = 0

Friends

class ParticleFilterPlacement


Constructor & Destructor Documentation

BLCMDparticlefilter::BLCMDparticlefilter (  ) 

Default constructor.

References BLCMDTYPE_ELEMENT, color, decay, height, innerRadius, keep, kill, length, material, maxStep, nWait, radius, referenceWait, BLCommand::registerCommand(), require, BLCommand::setDescription(), BLCommand::setSynopsis(), solid, steppingVerbose, and width.

Referenced by clone(), and command().

00167                                          : BLElement(), decaySet(), killSet(),
00168                                                                 keepSet()
00169 {
00170         registerCommand(BLCMDTYPE_ELEMENT);
00171         setSynopsis("Will kill particles from a list, or force particles to decay.");
00172         setDescription("A particlefilter will force the decay of certain particles\n"
00173                 "when they enter the physical volume of the element.\n"
00174                 "The list of affected particles is in the 'decay' argument,\n"
00175                 "and the normal Decay process is disabled for them.\n"
00176                 "In addition, the 'kill' argument is a list of particles to\n"
00177                 "kill when they enter the element, and the 'keep' argument\n"
00178                 "will kill all particles not named, if it is not empty.\n\n"
00179                 "require is an expression invloving track parameters that will "
00180                 "kill the track if it evaluates to zero (use a comparison "
00181                 "operator). The variables available are:\n"
00182                 "    x,y,z,t,Px,Py,Pz,Ptot,PDGid,TrkId,Parent\n"
00183                 "Units are mm, ns, MeV/c.\n\n"
00184                 "If nWait > 1, particles will not be killed until they hit "
00185                 "this element nwait times; this can be used to limit the "
00186                 "number of revolutions around a ring. Decays are unaffected "
00187                 "by nWait. referenceWait does the same for the reference "
00188                 "particle\n"
00189                 "The element can be either a cylinder or a box: set length and "
00190                 "radius, or set length and width and height.");
00191         
00192         // initial field values
00193         radius = -1.0;
00194         innerRadius = 0.0;
00195         height = 0.0;
00196         width = 0.0;
00197         length = 1.0*mm;
00198         material = "";
00199         color = "1,1,1";
00200         solid = 0;
00201         maxStep = -99.0;
00202         decay = "";
00203         kill = "";
00204         keep = "";
00205         nWait = 1;
00206         referenceWait = 1;
00207         require = "";
00208         steppingVerbose = 0;
00209 }

BLCMDparticlefilter::~BLCMDparticlefilter (  )  [virtual]

Destructor.

00212 {
00213 }

BLCMDparticlefilter::BLCMDparticlefilter ( const BLCMDparticlefilter r  ) 

Copy constructor.

References color, decay, height, innerRadius, keep, kill, length, material, maxStep, nWait, radius, referenceWait, require, solid, steppingVerbose, and width.

00215                                                                      : BLElement(r),
00216                 decaySet(r.decaySet), killSet(r.killSet), keepSet(r.keepSet)
00217 {
00218         radius = r.radius;
00219         innerRadius = r.innerRadius;
00220         height = r.height;
00221         width = r.width;
00222         length = r.length;
00223         material = r.material;
00224         color = r.color;
00225         solid = 0;
00226         maxStep = r.maxStep;
00227         decay = r.decay;
00228         kill = r.kill;
00229         keep = r.keep;
00230         nWait = r.nWait;
00231         referenceWait = r.referenceWait;
00232         require = r.require;
00233         steppingVerbose = r.steppingVerbose;
00234 }


Member Function Documentation

void BLCMDparticlefilter::setup (  )  [private]

References decay, decaySet, keep, keepSet, kill, killSet, and BLCommand::splitString().

Referenced by command().

00265 {
00266         // get set of particle names into decaySet
00267         std::vector<G4String> v=splitString(decay,',',true);
00268         for(unsigned i=0; i<v.size(); ++i) {
00269                 if(v[i].size() == 0) continue;
00270                 decaySet.insert(v[i]);
00271         }
00272 /***
00273         unsigned int place, next;
00274         for(place=next=0; next<decay.size(); place=next+1) {
00275                 next = decay.find(",",place);
00276                 G4String p;
00277                 if(next < decay.size())
00278                         p = decay.substr(place,next-place);
00279                 else
00280                         p = decay.substr(place);
00281                 if(p.size() == 0) break;
00282                 decaySet.insert(p);
00283         }
00284 ***/
00285         // get set of particle names into killSet
00286         v=splitString(kill,',',true);
00287         for(unsigned i=0; i<v.size(); ++i) {
00288                 if(v[i].size() == 0) continue;
00289                 killSet.insert(v[i]);
00290         }
00291 /***
00292         for(place=next=0; next<kill.size(); place=next+1) {
00293                 next = kill.find(",",place);
00294                 G4String p;
00295                 if(next < kill.size())
00296                         p = kill.substr(place,next-place);
00297                 else
00298                         p = kill.substr(place);
00299                 if(p.size() == 0) break;
00300                 killSet.insert(p);
00301         }
00302 ***/
00303         // get set of particle names into keepSet
00304         v=splitString(keep,',',true);
00305         for(unsigned i=0; i<v.size(); ++i) {
00306                 if(v[i].size() == 0) continue;
00307                 keepSet.insert(v[i]);
00308         }
00309 /***
00310         for(place=next=0; next<keep.size(); place=next+1) {
00311                 next = keep.find(",",place);
00312                 G4String p;
00313                 if(next < keep.size())
00314                         p = keep.substr(place,next-place);
00315                 else
00316                         p = keep.substr(place);
00317                 if(p.size() == 0) break;
00318                 keepSet.insert(p);
00319         }
00320 ***/
00321         // if keep is nonempty, add decaySet to it (they will decay)
00322         if(keepSet.size() > 0) {
00323                 std::set<G4String>::iterator i;
00324                 for(i=decaySet.begin(); i!=decaySet.end(); ++i)
00325                         keepSet.insert(*i);
00326         }
00327 }

BLElement* BLCMDparticlefilter::clone (  )  [inline, virtual]

clone()

Implements BLElement.

References BLCMDparticlefilter().

00085 { return new BLCMDparticlefilter(*this); }

G4String BLCMDparticlefilter::commandName (  )  [inline, virtual]

commandName() returns "particlefilter".

Implements BLCommand.

00091 { return "particlefilter"; }

int BLCMDparticlefilter::command ( BLArgumentVector argv,
BLArgumentMap namedArgs 
) [virtual]

command() implements the particlefilter command.

Implements BLCommand.

References BLCMDparticlefilter(), BLParam::getDouble(), BLParam::getInt(), BLCommand::getMaterial(), BLCommand::handleNamedArgs(), material, maxStep, Param, BLCommand::print(), BLCommand::printError(), BLElement::setName(), setup(), and steppingVerbose.

00237 {
00238         if(argv.size() != 1) {
00239                 printError("particlefilter: Invalid command, must have name");
00240                 return -1;
00241         }
00242 
00243         if(argv[0] == "default") {
00244                 return defaultParticleFilter.handleNamedArgs(namedArgs);
00245         }
00246 
00247         BLCMDparticlefilter *t = new BLCMDparticlefilter(defaultParticleFilter);
00248         t->setName(argv[0]);
00249         t->steppingVerbose = Param.getInt("steppingVerbose");
00250         int retval = t->handleNamedArgs(namedArgs);
00251 
00252         if(t->maxStep < 0.0) t->maxStep = Param.getDouble("maxStep");
00253 
00254         t->setup();
00255 
00256         // check material exists
00257         if(t->material.size() > 0) getMaterial(t->material,false);
00258 
00259         t->print(argv[0]);
00260 
00261         return retval;
00262 }

void BLCMDparticlefilter::defineNamedArgs (  )  [virtual]

defineNamedArgs() defines the named arguments for the command.

Reimplemented from BLCommand.

References BLCommand::argDouble(), BLCommand::argInt(), BLCommand::argString(), color, decay, height, innerRadius, keep, kill, length, material, maxStep, nWait, radius, referenceWait, require, steppingVerbose, and width.

00330 {
00331         argDouble(radius,"radius","The radius of the cylindrical particlefilter (mm).");
00332         argDouble(innerRadius,"innerRadius","The inner radius of the cylindrical particlefilter (0 mm, solid).");
00333         argDouble(height,"height","The height of the rectangular particlefilter (mm).");
00334         argDouble(width,"width","The width of the rectangular particlefilter (mm).");
00335         argDouble(length,"length","The length of the particlefilter (mm).");
00336         argDouble(maxStep,"maxStep","The maximum stepsize in the element (mm).");
00337         argString(material,"material","The material (default=parent's).");
00338         argString(color,"color","The color of the particlefilter (white).");
00339         argString(decay,"decay","A comma-separated list of particle names to decay.",false);
00340         argString(kill,"kill","A comma-separated list of particle names to kill.",false);
00341         argString(keep,"keep","A comma-separated list of particle names to keep.",false);
00342         argInt(nWait,"nWait","Intersection # to do the kill (default = 1)");
00343         argInt(referenceWait,"referenceWait","Intersection # for reference (default = 1)");
00344         argString(require,"require","Expression which will kill the track if zero.",false);
00345         argInt(steppingVerbose,"steppingVerbose","Nonzero to display track kills.");
00346         argString(decay,"decays","Synonym for decay.",false);
00347 
00348         if(radius > 0.0)
00349                 width = height = 2.0 * radius;
00350 }

void BLCMDparticlefilter::construct ( G4RotationMatrix *  relativeRotation,
G4ThreeVector  relativePosition,
G4LogicalVolume *  parent,
G4String  parentName,
G4RotationMatrix *  parentRotation,
G4ThreeVector  parentPosition 
) [virtual]

construct() will construct the particlefilter.

Implements BLElement.

References color, BLParam::getDouble(), BLCommand::getMaterial(), BLElement::getName(), BLManager::getObject(), BLCommand::getVisAttrib(), height, innerRadius, instance, keepSet, killSet, length, material, maxStep, minStep, Param, ParticleFilterPlacement, BLCommand::printError(), radius, BLManager::registerSteppingAction(), BLManager::registerTrackingAction(), require, solid, and width.

00358 {
00359         if(maxStep < 0.0) maxStep = Param.getDouble("maxStep");
00360 
00361         // ensure at least two steps will be taken inside the physical volume
00362         G4double minStep = (width<height ? width : height);
00363         minStep = (minStep<length ? minStep : length);
00364         minStep = minStep/2.0;
00365         maxStep = (minStep<maxStep ? minStep : maxStep);
00366 
00367         G4Material *mat;
00368         if(material != "")
00369                 mat = getMaterial(material);
00370         else
00371                 mat = parent->GetMaterial();
00372 
00373         G4String thisname = parentName+getName();
00374 
00375         if(!solid) {
00376                 if(radius > 0.0) {
00377                         solid = new G4Tubs(thisname+"Tubs", innerRadius, radius,
00378                                         length/2.0, 0.0, 2.0*pi);
00379                 } else if(height > 0.0 && width > 0.0) {
00380                         solid = new G4Box(thisname+"Box",width/2.0,
00381                                         height/2.0,length/2.0);
00382                 } else {
00383                         printError("particlefilter::construct %s INVALID - no "
00384                                 "radius or height&width",thisname.c_str());
00385                         return;
00386                 }
00387         }
00388         G4LogicalVolume *lv = new G4LogicalVolume(solid,mat, thisname+"LogVol");
00389         lv->SetVisAttributes(getVisAttrib(color));
00390         if(maxStep < 0.0) maxStep = Param.getDouble("maxStep");
00391         lv->SetUserLimits(new G4UserLimits(maxStep));
00392 
00393         // geant4 rotation convention is backwards from g4beamline
00394         G4RotationMatrix *g4rot = 0;
00395         if(relativeRotation)
00396                 g4rot = new G4RotationMatrix(relativeRotation->inverse());
00397 
00398         G4VPhysicalVolume *vol = new G4PVPlacement(g4rot, relativePosition,
00399                                         lv,thisname,parent,false,0);
00400 
00401         char processname[64];
00402         sprintf(processname,"BLCMDparticlefilter%d",++instance);
00403         ParticleFilterPlacement *place = 
00404                         new ParticleFilterPlacement(this,vol,processname);
00405 
00406         if(killSet.size() > 0 || keepSet.size() > 0 || require != "") {
00407             BLManager::getObject()->registerTrackingAction(place);
00408             BLManager::getObject()->registerSteppingAction(vol,place);
00409         }
00410 
00411         printf("BLCMDparticlefilter::Construct %s parent=%s\n",
00412                 thisname.c_str(),parentName.c_str());
00413 }

virtual G4double BLCMDparticlefilter::getLength (  )  [inline, virtual]

getLength() returns this element's Length along the Z axis.

Implements BLElement.

References length.

00108 { return length; }

virtual G4double BLCMDparticlefilter::getWidth (  )  [inline, virtual]

getWidth() returns this element's Width along the X axis.

Implements BLElement.

References width.

00111 { return width; }

virtual G4double BLCMDparticlefilter::getHeight (  )  [inline, virtual]

getHeight() returns this element's height along the Y axis.

Implements BLElement.

References height.

00114 { return height; }

virtual G4bool BLCMDparticlefilter::isOK (  )  [inline, virtual]

isOK() returns true.

Implements BLElement.

00117 { return true; }

void BLCMDparticlefilter::generatePoints ( int  npoints,
std::vector< G4ThreeVector > &  v 
) [virtual]

generatePoints() from BLElement

Implements BLElement.

References BLElement::generateBox(), BLElement::generateTubs(), height, innerRadius, length, radius, and width.

00416 {
00417         if(radius > 0.0)
00418                 generateTubs(npoints, innerRadius, radius, 0.0, 360.0*deg, length, v);
00419         else
00420                 generateBox(npoints,width,height,length,v);
00421 }

G4bool BLCMDparticlefilter::isOutside ( G4ThreeVector &  local,
G4double  tolerance 
) [virtual]

isOutside() from BLElement

Implements BLElement.

References height, innerRadius, length, radius, and width.

00424 {
00425         if(radius > 0.0) {
00426                 G4double r = sqrt(local[0]*local[0]+local[1]*local[1]);
00427                 return r > radius-tolerance || r < innerRadius+tolerance ||
00428                         fabs(local[2]) > length/2.0-tolerance;
00429         } else {
00430                 return fabs(local[0]) > width/2.0-tolerance ||
00431                         fabs(local[1]) > height/2.0-tolerance ||
00432                         fabs(local[2]) > length/2.0-tolerance;
00433         }
00434 }


Friends And Related Function Documentation

friend class ParticleFilterPlacement [friend]

Referenced by construct().


Member Data Documentation

G4double BLCMDparticlefilter::radius [private]

G4double BLCMDparticlefilter::height [private]

G4double BLCMDparticlefilter::width [private]

G4double BLCMDparticlefilter::length [private]

G4String BLCMDparticlefilter::material [private]

G4String BLCMDparticlefilter::color [private]

G4VSolid* BLCMDparticlefilter::solid [private]

Referenced by BLCMDparticlefilter(), and construct().

G4double BLCMDparticlefilter::maxStep [private]

G4String BLCMDparticlefilter::decay [private]

G4String BLCMDparticlefilter::kill [private]

G4String BLCMDparticlefilter::keep [private]

G4int BLCMDparticlefilter::nWait [private]

G4String BLCMDparticlefilter::require [private]

std::set<G4String> BLCMDparticlefilter::decaySet [private]

std::set<G4String> BLCMDparticlefilter::killSet [private]

std::set<G4String> BLCMDparticlefilter::keepSet [private]

int BLCMDparticlefilter::instance = 0 [static, private]

Referenced by construct().


The documentation for this class was generated from the following file:
g4beamline