TimeNTuple Class Reference

Inheritance diagram for TimeNTuple:

BLManager::SteppingAction

List of all members.


Detailed Description

class TimeNTuple implements an NTuple for a BLCMDtimentuple.

Public Member Functions

 TimeNTuple (G4String _name, G4double _time, G4String format, G4String filename, G4String _require, BLCoordinateType _coordinateType)
 constructor.
void UserSteppingAction (const G4Step *step)
 UserSteppingAction() from BLManager::SteppingAction.

Private Attributes

G4double time
BLNTuplentuple
G4String require
BLEvaluatoreval
BLCoordinateType coordinateType

Constructor & Destructor Documentation

TimeNTuple::TimeNTuple ( G4String  _name,
G4double  _time,
G4String  format,
G4String  filename,
G4String  _require,
BLCoordinateType  _coordinateType 
) [inline]

constructor.

References coordinateType, BLNTuple::create(), eval, ntuple, require, time, and TrackFields.

00087                                                           : 
00088                                                 BLManager::SteppingAction()
00089         { time = _time; 
00090           ntuple = BLNTuple::create(format,"TimeNTuple",_name,TrackFields,
00091                                                                 filename);
00092           require = _require;
00093           eval = 0;
00094           if(require != "") eval = new BLEvaluator();
00095           coordinateType = _coordinateType;
00096         }


Member Function Documentation

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

UserSteppingAction() from BLManager::SteppingAction.

Implements BLManager::SteppingAction.

References BLNTuple::appendRow(), BEAM, BLCOORD_REFERENCE, coordinateType, eval, BLEvaluator::evaluate(), BLCoordinates::getCoords(), BLManager::getExternalParentID(), BLManager::getExternalTrackID(), BLManager::getObject(), BLCoordinates::getRotation(), BLManager::getState(), BLEvaluator::isOK(), BLCoordinates::isValid(), NTrackFields, ntuple, BLCommand::printError(), require, BLEvaluator::setTrackVariables(), SPECIAL, and time.

00192 {
00193         if(BLManager::getObject()->getState() == SPECIAL) return;
00194 
00195         // only use reference coordinates when they are valid
00196         BLManagerState state = BLManager::getObject()->getState();
00197         if(coordinateType == BLCOORD_REFERENCE && state != BEAM) return;
00198 
00199         G4Track *track = step->GetTrack();
00200         if(eval) {
00201                 eval->setTrackVariables(track,coordinateType);
00202                 G4double v = eval->evaluate(require);
00203                 if(!eval->isOK()) {
00204                         BLCommand::printError("beamlossntuple: invalid expression require='%s'\n",
00205                                 require.c_str());
00206                         delete eval;
00207                         eval = 0;
00208                         v = 1.0;
00209                 }
00210                 if(v == 0.0) return;
00211         }
00212 
00213         // get basic physical-volume info
00214         G4StepPoint *prePoint = step->GetPreStepPoint();
00215         if(!prePoint) return;
00216         G4StepPoint *postPoint = step->GetPostStepPoint();
00217         if(!postPoint) return;
00218 
00219         // make a row only if the pre- and post-step times bracket desired time
00220         // and protect against ntuple==null
00221         G4double preTime = prePoint->GetGlobalTime();
00222         G4double postTime = postPoint->GetGlobalTime();
00223         if(preTime > time || postTime < time || !ntuple) return;
00224         
00225         G4RunManager* runmgr = G4RunManager::GetRunManager();
00226         const G4Event* event = runmgr->GetCurrentEvent();
00227         int evId = event->GetEventID();
00228         G4ThreeVector position = track->GetPosition();
00229         G4double time = track->GetGlobalTime();
00230         G4ThreeVector momentum = track->GetMomentum();
00231 
00232         // Linear interpolation
00233         G4double dt = step->GetDeltaTime();
00234         G4double err = this->time - time;
00235         if(dt > 0.001*ns  && fabs(err) <= dt*1.5) {
00236                 time += err;
00237                 position += (err/dt)*step->GetDeltaPosition();
00238                 momentum += (err/dt)*step->GetDeltaMomentum();
00239         }
00240 
00241         // transform to desired coordinates, if available
00242         BLCoordinates *c = (BLCoordinates *)track->GetUserInformation();
00243         if(c && c->isValid()) {
00244                 c->getCoords(coordinateType,position);
00245                 momentum = c->getRotation() * momentum;
00246         }
00247 
00248         float data[NTrackFields];
00249         data[0] = position[0]/mm;               // x (mm)
00250         data[1] = position[1]/mm;               // y (mm)
00251         data[2] = position[2]/mm;               // z (mm)
00252         data[3] = momentum[0]/MeV;              // Px (MeV/c)
00253         data[4] = momentum[1]/MeV;              // Py (MeV/c)
00254         data[5] = momentum[2]/MeV;              // Pz (MeV/c)
00255         data[6] = time/ns;                      // t (ns)
00256         data[7] = track->GetDefinition()->GetPDGEncoding();
00257         data[8] = evId;                         // Event ID
00258         data[9] = BLManager::getObject()->getExternalTrackID(track);
00259         data[10] = BLManager::getObject()->getExternalParentID(track);
00260         data[11] = track->GetWeight();          // Weight
00261 
00262         ntuple->appendRow(data,NTrackFields);
00263 }


Member Data Documentation

G4double TimeNTuple::time [private]

Referenced by TimeNTuple(), and UserSteppingAction().

Referenced by TimeNTuple(), and UserSteppingAction().

G4String TimeNTuple::require [private]

Referenced by TimeNTuple(), and UserSteppingAction().

Referenced by TimeNTuple(), and UserSteppingAction().

Referenced by TimeNTuple(), and UserSteppingAction().


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