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.

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


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