TraceNTuple Class Reference

List of all members.


Detailed Description

class TraceNTuple implements an NTuple for tracing a track.

Public Member Functions

 TraceNTuple (G4String format, G4String name, G4String filename, BLCoordinateType _coordinateType)
 constructor.
virtual ~TraceNTuple ()
 destructor.
void appendTrace (G4Track *track)
 appendTrace() will append a track to the Trace NTuple.
void close ()
 close() will close the NTuple.
void annotate (G4String line)
 annotate() will add an annotation to ASCII NTuples.

Private Attributes

BLNTuplentuple
BLCoordinateType coordinateType
BLManagermanager
BLRunManagerrunManager

Constructor & Destructor Documentation

TraceNTuple::TraceNTuple ( G4String  format,
G4String  name,
G4String  filename,
BLCoordinateType  _coordinateType 
)

constructor.

References coordinateType, BLNTuple::create(), BLRunManager::getObject(), BLManager::getObject(), manager, ntuple, runManager, and TraceFields.

00361 {
00362         ntuple = BLNTuple::create(format,"Trace",name,TraceFields,filename);
00363         coordinateType = _coordinateType;
00364         manager = BLManager::getObject();
00365         runManager = BLRunManager::getObject();
00366 }

TraceNTuple::~TraceNTuple (  )  [virtual]

destructor.

References ntuple.

00369 {
00370         if(ntuple) delete ntuple;
00371         ntuple = 0;
00372 }


Member Function Documentation

void TraceNTuple::appendTrace ( G4Track *  track  ) 

appendTrace() will append a track to the Trace NTuple.

References BLNTuple::appendRow(), B, BEAM, BLCOORD_CENTERLINE, coordinateType, E, BLCoordinates::getCoords(), BLManager::getExternalParentID(), BLManager::getExternalTrackID(), BLGlobalField::GetFieldValue(), BLGlobalField::getObject(), BLCoordinates::getRotation(), BLManager::getState(), BLCoordinates::isValid(), manager, NTraceFields, ntuple, and runManager.

Referenced by BLCMDtrace::collectiveStep(), and BLCMDtrace::UserSteppingAction().

00375 {
00376         const G4Event* event = runManager->GetCurrentEvent();
00377         int evNum = event->GetEventID();
00378         G4ThreeVector position = track->GetPosition();
00379         G4double time = track->GetGlobalTime();
00380         G4ThreeVector momentum = track->GetMomentum();
00381 
00382         G4double point[4], field[6];
00383         point[0] = position[0];
00384         point[1] = position[1];
00385         point[2] = position[2];
00386         point[3] = time;
00387         BLGlobalField::getObject()->GetFieldValue(point,field);
00388         G4ThreeVector B(field[0],field[1],field[2]);
00389         G4ThreeVector E(field[3],field[4],field[5]);
00390 
00391         // transform to desired coordinates, if appropriate
00392         BLCoordinates *c = (BLCoordinates *)track->GetUserInformation();
00393         if(c && c->isValid()) {
00394                 BLManagerState state = manager->getState();
00395                 c->getCoords((state==BEAM?coordinateType:BLCOORD_CENTERLINE),
00396                                                                 position);
00397                 momentum = c->getRotation() * momentum;
00398                 B = c->getRotation() * B;
00399                 E = c->getRotation() * E;
00400         }
00401 
00402         float data[NTraceFields];
00403         data[0] = position[0]/mm;               // x (mm)
00404         data[1] = position[1]/mm;               // y (mm)
00405         data[2] = position[2]/mm;               // z (mm)
00406         data[3] = momentum[0]/MeV;              // Px (MeV/c)
00407         data[4] = momentum[1]/MeV;              // Py (MeV/c)
00408         data[5] = momentum[2]/MeV;              // Pz (MeV/c)
00409         data[6] = time/ns;                      // t (ns)
00410         data[7] = track->GetDefinition()->GetPDGEncoding();
00411         data[8] = evNum;                                // Event ID
00412         data[9] = manager->getExternalTrackID(track);
00413         data[10] = manager->getExternalParentID(track);
00414         data[11] = track->GetWeight();          // Weight
00415         data[12] = B[0]/tesla;                  // Bx
00416         data[13] = B[1]/tesla;                  // By
00417         data[14] = B[2]/tesla;                  // Bz
00418         data[15] = E[0]/(megavolt/meter);       // Ex
00419         data[16] = E[1]/(megavolt/meter);       // Ey
00420         data[17] = E[2]/(megavolt/meter);       // Ez
00421 
00422         ntuple->appendRow(data,NTraceFields);
00423 }

void TraceNTuple::close (  )  [inline]

void TraceNTuple::annotate ( G4String  line  )  [inline]

annotate() will add an annotation to ASCII NTuples.

References BLNTuple::annotate(), and ntuple.

Referenced by BLCMDtrace::PostUserTrackingAction(), and BLCMDtrace::PreUserTrackingAction().

00065 { if(ntuple) ntuple->annotate(line); }


Member Data Documentation

Referenced by appendTrace(), and TraceNTuple().

Referenced by appendTrace(), and TraceNTuple().

Referenced by appendTrace(), and TraceNTuple().


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