BLVirtualDetectorNTuple Class Reference

Inheritance diagram for BLVirtualDetectorNTuple:

BLManager::SteppingAction

List of all members.


Detailed Description

class BLVirtualDetectorNTuple implements an NTuple for a BLCMDvirtualdetector.

Public Member Functions

 BLVirtualDetectorNTuple (G4String type, G4String category, G4String _name, G4VPhysicalVolume *pv, int noSingles, G4String filename, G4String _require, BLCoordinateType _coordinateType, int _kill)
 constructor.
bool needsReference ()
void UserSteppingAction (const G4Step *step)
 UserSteppingAction() from BLManager::SteppingAction.

Private Attributes

G4String name
G4VPhysicalVolume * thisVol
BLCoordinateType coordinateType
BLNTuplentuple
int nNTuples
G4String require
BLEvaluatoreval
int kill
bool verbose

Friends

class BLCMDvirtualdetector

Constructor & Destructor Documentation

BLVirtualDetectorNTuple::BLVirtualDetectorNTuple ( G4String  type,
G4String  category,
G4String  _name,
G4VPhysicalVolume *  pv,
int  noSingles,
G4String  filename,
G4String  _require,
BLCoordinateType  _coordinateType,
int  _kill 
)

constructor.

References coordinateType, BLNTuple::create(), eval, BLManager::getObject(), BLManager::getSteppingVerbose(), kill, name, nNTuples, ntuple, require, BLCMDntuple::setDetector(), thisVol, TrackFields, and verbose.

00374 {
00375         name = _name;
00376         thisVol = pv;
00377         if(noSingles == 0)
00378                 ntuple = BLNTuple::create(type,category,name,TrackFields,
00379                                                                 filename);
00380         else
00381                 ntuple = 0;
00382         coordinateType = _coordinateType;
00383         require = _require;
00384         eval = 0;
00385         if(require != "") eval = new BLEvaluator();
00386         kill = _kill;
00387         verbose = BLManager::getObject()->getSteppingVerbose() > 0;
00388 
00389         // define this virtualdetector to the ntuple command
00390         nNTuples = BLCMDntuple::setDetector(this,name,TrackFields);
00391 }


Member Function Documentation

bool BLVirtualDetectorNTuple::needsReference (  )  [inline]

References BLNTuple::needsReference(), and ntuple.

Referenced by BLCMDvirtualdetector::construct().

00161                 { return ntuple ? ntuple->needsReference() : false; }

void BLVirtualDetectorNTuple::UserSteppingAction ( const G4Step *  step  )  [virtual]

UserSteppingAction() from BLManager::SteppingAction.

Implements BLManager::SteppingAction.

References BLCMDntuple::appendRow(), BLNTuple::appendRow(), B, BEAM, BLCOORD_REFERENCE, coordinateType, E, eval, BLEvaluator::evaluate(), BLCoordinates::getCoords(), BLManager::getExternalParentID(), BLManager::getExternalTrackID(), BLGlobalField::GetFieldValue(), BLNTuple::getNData(), BLGlobalField::getObject(), BLManager::getObject(), BLCoordinates::getRotation(), BLManager::getState(), BLEvaluator::isOK(), BLCoordinates::isValid(), kill, name, nNTuples, NTrackFields, ntuple, BLCommand::printError(), require, BLEvaluator::setTrackVariables(), SPECIAL, thisVol, and verbose.

00394 {
00395         if(BLManager::getObject()->getState() == SPECIAL) return;
00396 
00397         // only use reference coordinates when they are valid
00398         BLManagerState state = BLManager::getObject()->getState();
00399         if(coordinateType == BLCOORD_REFERENCE && state != BEAM) return;
00400 
00401         G4Track *track = step->GetTrack();
00402         if(eval) {
00403                 eval->setTrackVariables(track,coordinateType);
00404                 G4double v = eval->evaluate(require);
00405                 if(!eval->isOK()) {
00406                         BLCommand::printError("beamlossntuple: invalid expression require='%s'\n",
00407                                 require.c_str());
00408                         delete eval;
00409                         eval = 0;
00410                         v = 1.0;
00411                 }
00412                 if(v == 0.0) return;
00413         }
00414 
00415         // get basic physical-volume info
00416         G4StepPoint *prePoint = step->GetPreStepPoint();
00417         if(!prePoint) return;
00418         G4VPhysicalVolume *preVol = prePoint->GetPhysicalVolume();
00419         if(!preVol) return;
00420         G4StepPoint *postPoint = step->GetPostStepPoint();
00421         if(!postPoint) return;
00422         G4VPhysicalVolume *postVol = postPoint->GetPhysicalVolume();
00423         if(!postVol) return;
00424         
00425         // return if not entering thisVol
00426         if(preVol == postVol || postVol != thisVol) return;
00427 
00428         G4RunManager* runmgr = G4RunManager::GetRunManager();
00429         const G4Event* event = runmgr->GetCurrentEvent();
00430         int evId = event->GetEventID();
00431         G4ThreeVector position = track->GetPosition();
00432         G4double time = track->GetGlobalTime();
00433         G4ThreeVector momentum = track->GetMomentum();
00434 
00435         // get # data fields, and get B and E if needed
00436         unsigned ndata = (ntuple ? ntuple->getNData() : NTrackFields);
00437         G4double field[6];
00438         if(ndata > NTrackFields) {
00439                 G4double point[4];
00440                 point[0] = position[0];
00441                 point[1] = position[1];
00442                 point[2] = position[2];
00443                 point[3] = time;
00444                 BLGlobalField::getObject()->GetFieldValue(point,field);
00445         }
00446 
00447         // transform to desired coordinates, if available
00448         BLCoordinates *c = (BLCoordinates *)track->GetUserInformation();
00449         if(c && c->isValid()) {
00450                 c->getCoords(coordinateType,position);
00451                 momentum = c->getRotation() * momentum;
00452         } else {
00453                 printf("BLCMDvirtualdetector::SteppingAction: track has no "
00454                         "BLCoordinates object\n");
00455         }
00456 
00457         float data[18];
00458         assert(ndata <= sizeof(data)/sizeof(data[0]));
00459         data[0] = position[0]/mm;               // x (mm)
00460         data[1] = position[1]/mm;               // y (mm)
00461         data[2] = position[2]/mm;               // z (mm)
00462         data[3] = momentum[0]/MeV;              // Px (MeV/c)
00463         data[4] = momentum[1]/MeV;              // Py (MeV/c)
00464         data[5] = momentum[2]/MeV;              // Pz (MeV/c)
00465         data[6] = time/ns;                      // t (ns)
00466         data[7] = track->GetDefinition()->GetPDGEncoding();
00467         data[8] = evId;                         // Event ID
00468         data[9] = BLManager::getObject()->getExternalTrackID(track);
00469         data[10] = BLManager::getObject()->getExternalParentID(track);
00470         data[11] = track->GetWeight();          // Weight
00471         if(ndata > NTrackFields) {
00472                 G4ThreeVector B(field[0],field[1],field[2]);
00473                 G4ThreeVector E(field[3],field[4],field[5]);
00474                 if(c && c->isValid()) {
00475                         B = c->getRotation() * B;
00476                         E = c->getRotation() * E;
00477                 }
00478                 data[12] = B[0]/tesla;
00479                 data[13] = B[1]/tesla;
00480                 data[14] = B[2]/tesla;
00481                 data[15] = E[0]/(megavolt/meter);
00482                 data[16] = E[1]/(megavolt/meter);
00483                 data[17] = E[2]/(megavolt/meter);
00484         }
00485 
00486         if(ntuple) ntuple->appendRow(data,ndata);
00487 
00488         // insert data into BLCMDntuple-s
00489         if(nNTuples > 0) BLCMDntuple::appendRow(this,data,NTrackFields);
00490 
00491         if(kill) {
00492                 track->SetTrackStatus(fStopAndKill);
00493                 if(verbose) printf("Track killed by '%s' with kill=1\n",
00494                         name.c_str());
00495         }
00496 }


Friends And Related Function Documentation

friend class BLCMDvirtualdetector [friend]


Member Data Documentation

G4String BLVirtualDetectorNTuple::name [private]

G4VPhysicalVolume* BLVirtualDetectorNTuple::thisVol [private]


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