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
int noSingles
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, noSingles, ntuple, require, thisVol, TrackFields, and verbose.

00378 {
00379         name = _name;
00380         thisVol = pv;
00381         noSingles = _noSingles;
00382         ntuple = BLNTuple::create(type,category,name,TrackFields,filename);
00383         coordinateType = _coordinateType;
00384         require = _require;
00385         eval = 0;
00386         if(require != "") eval = new BLEvaluator();
00387         kill = _kill;
00388         verbose = BLManager::getObject()->getSteppingVerbose() > 0;
00389 }


Member Function Documentation

bool BLVirtualDetectorNTuple::needsReference (  )  [inline]

References BLNTuple::needsReference(), and ntuple.

Referenced by BLCMDvirtualdetector::construct().

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

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

UserSteppingAction() from BLManager::SteppingAction.

Implements BLManager::SteppingAction.

References BLNTuple::appendRow(), B, BEAM, BLCOORD_REFERENCE, coordinateType, BLNTuple::doCallbacks(), 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, noSingles, NTrackFields, ntuple, BLCommand::printError(), require, BLEvaluator::setTrackVariables(), SPECIAL, thisVol, and verbose.

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


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