Skip to content

Commit c58ee23

Browse files
bazinskibazinski
committed
mid way through building tree in the raw data manager
1 parent c4f3afa commit c58ee23

4 files changed

Lines changed: 132 additions & 39 deletions

File tree

‎Detectors/TRD/qc/include/TRDQC/CoordinateTransformer.h‎

Lines changed: 9 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -151,7 +151,9 @@ class TrackSegment
151151
void setCollisionId(int id) { mCollisionId=id;}
152152

153153
int getTriggerTime() const { return mTriggerTime;}
154-
void setTriggerTime(float ttime) { mTriggerTime=ttime;}
154+
void setTriggerTime(float time) { mTriggerTime=time;}
155+
int getTrackTime() const { return mTrackTime;}
156+
void setTrackTime(float time) { mTrackTime=time;}
155157

156158
ChamberSpacePoint getStartPoint()const { return mStartPoint; }
157159
ChamberSpacePoint getEndPoint() const { return mEndPoint; }
@@ -177,6 +179,8 @@ class TrackSegment
177179
void setStartX(float x){ mStartX=x;}
178180
void setStartY(float y){ mStartY=y;}
179181
void setStartZ(float z){ mStartZ=z;}
182+
int getLayer(){ return mLayer;}
183+
void setLayer(int layer){mLayer=layer;}
180184

181185
protected:
182186
ChamberSpacePoint mStartPoint, mEndPoint;
@@ -186,6 +190,7 @@ class TrackSegment
186190
int mRefITSId;
187191
int mCollisionId{0};
188192
float mTriggerTime{0.0};
193+
float mTrackTime{0.0};
189194
float mSagita; // this is findable from the TrackID
190195
int mDetector;
191196
int mDetector2; // for those instances where the track finishes a layer in a different detector to that which is starts in.
@@ -196,7 +201,7 @@ class TrackSegment
196201
float mStartY{0.0};
197202
float mStartZ{0.0};
198203
float mAlpha{0.0};
199-
204+
int mLayer;
200205

201206

202207
};
@@ -243,6 +248,8 @@ class CoordinateTransformer
243248

244249
// recalculate given a new t0, cdrift and exb
245250
std::array<float, 3> RecalculateRCT(int det, float x, float y, float z, float t0, float vdrift, float exb);
251+
std::array<float, 3> RecalculateRCT(int det, o2::trd::ChamberSpacePoint& point, float t0, float vdrift, float exb);
252+
std::array<float, 3> RecalculateRCT(int det,o2::trd::ChamberSpacePoint& point);
246253
std::array<float, 3> RecalculateRCT(int det, float x, float y, float z);
247254

248255
/// Legacy, less accurate method to convert local spatial to row/column/time coordinate.

‎Detectors/TRD/qc/include/TRDQC/RawDataManager.h‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -242,6 +242,9 @@ class RawDataManager
242242

243243
float mAverageRadii[6] = {300.2f, 312.8f, 325.4f, 338.0f, 350.6f, 363.2f}; // used as default value in case no transformation matrix can be obtained
244244
//GPUTRDGeometry mGeogpu; // TRD geometry
245+
// only used for those times we write the track to tracklet matching to a root tree.
246+
TFile *mfile=nullptr;
247+
TTree *moutputtree=nullptr;
245248
//
246249
//
247250
//const int kNStacks=5;

‎Detectors/TRD/qc/src/CoordinateTransformer.cxx‎

Lines changed: 20 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -145,6 +145,26 @@ std::array<float, 3> CoordinateTransformer::RecalculateRCT(int det, float x, flo
145145
mVdrift=vdrift;
146146
return Local2RCT(det,x,y,z);
147147
}
148+
std::array<float, 3> CoordinateTransformer::RecalculateRCT(int det, ChamberSpacePoint& p, float t0, float vdrift, float exb)
149+
{
150+
float x,y,z;
151+
x=p.getX();
152+
y=p.getY();
153+
z=p.getZ();
154+
mT0=t0;
155+
mExB=exb;
156+
mVdrift=vdrift;
157+
return Local2RCT(det,x,y,z);
158+
}
159+
160+
std::array<float, 3> CoordinateTransformer::RecalculateRCT(int det, ChamberSpacePoint& p)
161+
{
162+
float x,y,z;
163+
x=p.getX();
164+
y=p.getY();
165+
z=p.getZ();
166+
return Local2RCT(det,x,y,z);
167+
}
148168

149169
std::array<float, 3> CoordinateTransformer::RecalculateRCT(int det, float x, float y, float z)
150170
{

‎Detectors/TRD/qc/src/RawDataManager.cxx‎

Lines changed: 100 additions & 37 deletions
Original file line numberDiff line numberDiff line change
@@ -369,7 +369,7 @@ std::vector<TrackSegment> RawDataSpan::makeMCTrackSegments()
369369
}
370370

371371
/// The RawDataManager constructor: connects all data files and sets up trees, readers etc.
372-
RawDataManager::RawDataManager(std::filesystem::path dir)
372+
RawDataManager::RawDataManager(std::filesystem::path dir, std::string treefname)
373373
{
374374

375375
if (!std::filesystem::exists(dir) || !std::filesystem::is_directory(dir)) {
@@ -555,8 +555,9 @@ RawDataManager::RawDataManager(std::filesystem::path dir)
555555
}
556556
mTransformer.init();
557557
mTransformer.setCalVdriftExB(calvdriftexb);
558-
559-
558+
LOGP(info,"Setup root tree for output");
559+
mfile = new TFile(treefname.c_str(),"RECREATE");
560+
moutputtree = new TTree("t","tracklets and tracksegments");
560561

561562
}
562563
/***********************************************************************************************/
@@ -932,6 +933,27 @@ bool RawDataManager::getYZAt(float xk, float b, float& y, float& z, o2::dataform
932933
return true;
933934
}
934935

936+
int RawDataManager::findNearestTracklet(o2::trd::TrackSegment& tracksegment)
937+
{
938+
auto triggertime = tracksegment.getTriggerTime();
939+
//loop through tracklets and find closest.
940+
//
941+
//
942+
//
943+
////// WE ARE HERE !!!!!!!
944+
for(auto& trgrec : mTriggerRecord){
945+
auto triggertime = getTriggerTime(trgrec, *mTFIDs, mTimeFrameNo - 1);
946+
if (!trackMatchesCollision(triggertime, tracksegment.getTrackTime())) {
947+
continue;
948+
}
949+
for(int trklt=trgrec.getFirstTracklet(); trklt<trgrec.getFirstTracklet()+trgrec.getNumberOfTracklets();++trklt){
950+
951+
}
952+
953+
954+
}
955+
956+
}
935957
int RawDataManager::propagateTrack(o2::dataformats::TrackTPCITS& track, float e, float maxStep, float triggertime, int& glbTrkltIdxOffset, int collisionId)
936958
{
937959
if (debugprint)
@@ -941,9 +963,33 @@ int RawDataManager::propagateTrack(o2::dataformats::TrackTPCITS& track, float e,
941963
int layerCount = 0;
942964
const int32_t nMaxChambersToSearch = 18;
943965
int32_t trkltIdxOffset = collisionId * (o2::trd::constants::NCHAMBER + 1); // offset for accessing mTrackletIndexArray for given collision
966+
std::vector<o2::trd::TrackSegment> btracksegments;
967+
std::vector<o2::trd::Tracklet64> btracklets;
968+
int timeframe;
969+
int eventno;
970+
int matched;
971+
float trackalpha,tracksnp,pT,trackxstart;
972+
float tracktime,trdtriggertime;
973+
int itsindex,tpcindex;
974+
int layers;
975+
moutputtree->Branch("segment", &btracksegments);
976+
moutputtree->Branch("tracklet", &btracklets);
977+
moutputtree->Branch("timeframe",&timeframe);
978+
//moutputtree->Branch("trackxstart",&trackxstart);
979+
moutputtree->Branch("alpha",&trackalpha);
980+
moutputtree->Branch("snp",&tracksnp);
981+
moutputtree->Branch("pt",&pT);
982+
moutputtree->Branch("itsindex",&itsindex);
983+
moutputtree->Branch("tpcindex",&tpcindex);
984+
moutputtree->Branch("event",&eventno);
985+
moutputtree->Branch("matched",&matched);
986+
moutputtree->Branch("triggertime",&trdtriggertime);
987+
moutputtree->Branch("tracktime",&tracktime);
988+
moutputtree->Branch("layers",&layers);
989+
//
944990
for (int32_t iLayer = 0; iLayer < 6; ++iLayer) {
945991
// nCurrHypothesis = 0;
946-
// if (debugprint)
992+
if (debugprint)
947993
LOGP(info, " $$$$ Layer : {} layercount: {}", iLayer,layerCount);
948994
const o2::trd::PadPlane* pad = mGeo->getPadPlane(iLayer, 0);
949995
float tilt = std::tan(std::numbers::pi / 180.f * pad->getTiltingAngle());
@@ -973,8 +1019,7 @@ int RawDataManager::propagateTrack(o2::dataformats::TrackTPCITS& track, float e,
9731019
float zShiftTrk = (track.getTimeMUS().getTimeStamp() - triggertime) * mTPCVdrift * side;
9741020

9751021
if (!isGeoFindable(track, iLayer, track.getAlpha(), zShiftTrk)) {
976-
LOGP(info,"Track not geofindable");
977-
1022+
// LOGP(info,"Track not geofindable");
9781023
continue;
9791024
}
9801025
layerCount++;
@@ -1006,7 +1051,6 @@ int RawDataManager::propagateTrack(o2::dataformats::TrackTPCITS& track, float e,
10061051
// LOGP(info,"chambers to search : {}",s);
10071052
// LOGP(info,"chambers to search : {}",aa);
10081053

1009-
// look for tracklets in chamber(s)
10101054
for (int32_t iDet = 0; iDet < nMaxChambersToSearch; iDet++) {
10111055
int32_t currDet = det[iDet];
10121056
if (currDet == -1) {
@@ -1023,11 +1067,11 @@ int RawDataManager::propagateTrack(o2::dataformats::TrackTPCITS& track, float e,
10231067
}
10241068
}
10251069
if (currSec != getSector(track.getAlpha())) {
1026-
// if (debugprint)
1070+
if (debugprint)
10271071
LOGP(info, "Track is in sector {} and we are in sector {}", getSector(track.getAlpha()), currSec);
10281072
continue;
10291073
}
1030-
// propagate track to radius of chamber
1074+
// propagate track to start radius of chamber drift start
10311075
const PadPlane* pp = mGeo->getPadPlane(currDet);
10321076
if (propagateToLayerX(track, mRdriftstart[currDet], 0.8f, 0.2f)) { // prop.propagateToX(mR[currDet], .8f, .2f)) {
10331077
// we are at the start of a layer now propagate to the outter radius and build a tracksegment for the voxel of an mcm. TODO voxel of a padrow
@@ -1078,6 +1122,7 @@ int RawDataManager::propagateTrack(o2::dataformats::TrackTPCITS& track, float e,
10781122
tracksegment.setEndPoint(ae);
10791123
tracksegment.setCollisionId(collisionId);
10801124
tracksegment.setTriggerTime(triggertime);
1125+
tracksegment.setTrackTime(track.getTimeMUS().getTimeStamp());
10811126
tracksegment.setRefTPCId(track.getRefTPC());
10821127
tracksegment.setRefITSId(track.getRefITS());
10831128
tracksegment.setPhi(track.getPhi());
@@ -1095,8 +1140,12 @@ int RawDataManager::propagateTrack(o2::dataformats::TrackTPCITS& track, float e,
10951140
if (debugprint)
10961141
LOGP(info, "TrackSegmentG x:y:z {:.2f}:{:.2f}:{:.2f} --> {:.2f}:{:.2f}:{:.2f}", localpointG.X(), localpointG.Y(), localpointG.Z(), localpointendG.X(), localpointendG.Y(), localpointendG.Z());
10971142
mITSTPCTracks_segments.push_back(tracksegment);
1143+
btracksegments.push_back(tracksegment);
1144+
//find nearest tracklet to this track segment.
1145+
int trackletindex=findNearestTracklet(tracksegment);
1146+
btracklet.push_back(mTracklets[trackletindex]);
10981147
} else {
1099-
LOGP(info, "Track could not be propagated to radius of chamber {}", currDet);
1148+
LOGP(info, "Track could not be propagated to radius of chamber {} which is layer: {}", currDet, iLayer);
11001149
}
11011150
} // chamber loop
11021151

@@ -1105,9 +1154,21 @@ int RawDataManager::propagateTrack(o2::dataformats::TrackTPCITS& track, float e,
11051154
} // end of layer loop
11061155
// if(layerCount>0) LOGP(info,"Layer count : {}",layerCount);
11071156
/****************************************************************************************/
1157+
//now write the tree event ..
1158+
timeframe=getTimeFrameNumber();
1159+
eventno=getEventNumber();
1160+
trackalpha=track.getAlpha();
1161+
tracksnp=track.getSnp();
1162+
pT=track.getPt();
1163+
itsindex=track.getRefITS();
1164+
tpcindex=track.getRefTPC();
1165+
tracktime=track.getTimeMUS().getTimeStamp();
1166+
trdtriggertime=triggertime;;
1167+
layers=layerCount;
1168+
11081169
// if (layerCount > 2)
11091170
// return true;
1110-
LOGP(info,"propagateTrack returning with layercount of {}",layerCount);
1171+
//LOGP(info,"propagateTrack returning with layercount of {}",layerCount);
11111172
return layerCount;
11121173
}
11131174
/***********************************************************************************************/
@@ -1120,6 +1181,7 @@ bool RawDataManager::buildTrackSegments(bool onlydigits)
11201181
//LOGP(info,"{} {}",__func__,__LINE__);
11211182
bool timeframehasdigits = false;
11221183
bool timeframehastracks = false;
1184+
11231185
//clear the tracksegments
11241186
mITSTPCTracks_segments.clear();
11251187
if (onlydigits) {
@@ -1198,36 +1260,37 @@ bool RawDataManager::buildTrackSegments(bool onlydigits)
11981260
auto ttrack = track;
11991261
trackcounter++;
12001262
if (trackMatchesCollision(triggertime, ttrack.getTimeMUS().getTimeStamp())) {
1201-
trackcounter++;
1202-
int numberOfLayers=0;
1203-
auto difftime = ttrack.getTimeMUS().getTimeStamp() - triggertime;
1204-
auto pt = track.getPt();
1205-
float timeWindow = 4.0; // time is within 20us
1206-
if (difftime < timeWindow && pt>0.5) {
1207-
// TracksForThisEvent.push_back(track);
1208-
// Track_pad_row_timebin.push_back(GeneratePadRowTimeBin(track));
1209-
numberOfTimeTrackMatched++;
1210-
tracktimewindowcounter++;
1211-
// if(pt>1.0){
1212-
// LOGP(info,"$$$$ Propagating track ..... {} {} track.x={} trigtime:{} tracktime:{} its:{} tpc:{} trackcounter:{} collionsId:{} timeframe:{} eventno:{}",__func__,__LINE__,ttrack.getX(),triggertime,ttrack.getTimeMUS().getTimeStamp(),
1213-
// (int)ttrack.getRefITS(),(int)ttrack.getRefTPC(),trackcounter,collisionId,mTimeFrameNo,mEventNo);
1214-
numberOfLayers=propagateTrack(ttrack, .8f, 2.f, triggertime, trackletstart, collisionId);
1215-
totalNumberOfLayers+=numberOfLayers;
1216-
if (numberOfLayers > -1) {
1217-
// LOGP(info," track propagated..... collid:{} timeframe:{} eventno:{}",collisionId,mTimeFrameNo,mEventNo);
1218-
goodtrackcounter++; // we have at least a singular layer 0,1,2,3,4,5;
1219-
} else {
1220-
badtrackcounter++;
1221-
// LOGP(info," track failed to propagate ..... collid:{} timeframe:{} eventno:{}",collisionId,mTimeFrameNo,mEventNo);
1222-
}
1223-
LOGP(info,"$$$$ Finished Propagating track ..... with track.x={} collid:{} timeframe:{} eventno:{} layerscount:{} layersfromprop:{}, trackpt:{}",ttrack.getX(),collisionId,mTimeFrameNo,mEventNo,totalNumberOfLayers, numberOfLayers, ttrack.getPt());
1224-
//}
1263+
continue;
1264+
}
1265+
trackcounter++;
1266+
int numberOfLayers=0;
1267+
auto difftime = ttrack.getTimeMUS().getTimeStamp() - triggertime;
1268+
auto pt = track.getPt();
1269+
float timeWindow = 4.0; // time is within 20us
1270+
if (difftime < timeWindow && pt>2.0) {
1271+
// TracksForThisEvent.push_back(track);
1272+
// Track_pad_row_timebin.push_back(GeneratePadRowTimeBin(track));
1273+
numberOfTimeTrackMatched++;
1274+
tracktimewindowcounter++;
1275+
// if(pt>1.0){
1276+
// LOGP(info,"$$$$ Propagating track ..... {} {} track.x={} trigtime:{} tracktime:{} its:{} tpc:{} trackcounter:{} collionsId:{} timeframe:{} eventno:{}",__func__,__LINE__,ttrack.getX(),triggertime,ttrack.getTimeMUS().getTimeStamp(),
1277+
// (int)ttrack.getRefITS(),(int)ttrack.getRefTPC(),trackcounter,collisionId,mTimeFrameNo,mEventNo);
1278+
numberOfLayers=propagateTrack(ttrack, .8f, 2.f, triggertime, trackletstart, collisionId);
1279+
totalNumberOfLayers+=numberOfLayers;
1280+
if (numberOfLayers > -1) {
1281+
// LOGP(info," track propagated..... collid:{} timeframe:{} eventno:{}",collisionId,mTimeFrameNo,mEventNo);
1282+
goodtrackcounter++; // we have at least a singular layer 0,1,2,3,4,5;
1283+
} else {
1284+
badtrackcounter++;
1285+
// LOGP(info," track failed to propagate ..... collid:{} timeframe:{} eventno:{}",collisionId,mTimeFrameNo,mEventNo);
12251286
}
1287+
//LOGP(info,"$$$$ Finished Propagating track ..... with track.x={} collid:{} timeframe:{} eventno:{} layerscount:{} layersfromprop:{}, trackpt:{}",ttrack.getX(),collisionId,mTimeFrameNo,mEventNo,totalNumberOfLayers, numberOfLayers, ttrack.getPt());
1288+
//}
12261289
}
12271290

12281291
collisionId++;
12291292
}
1230-
LOGP(info,"Tracks in collisionid {} this event {} out of a total of {} good tracks from {} tracks, with {}% good events, track matched {} times numberoflayers:{}", collisionId, goodtrackcounter, trackcounter, (*mITSTPCTracks).size(), (float)goodtrackcounter/(float)trackcounter*100,numberOfTimeTrackMatched,totalNumberOfLayers);
1293+
//LOGP(info,"Tracks in collisionid {} this event {} out of a total of {} good tracks from {} tracks, with {}% good events, track matched {} times numberoflayers:{}", collisionId, goodtrackcounter, trackcounter, (*mITSTPCTracks).size(), (float)goodtrackcounter/(float)trackcounter*100,numberOfTimeTrackMatched,totalNumberOfLayers);
12311294
/********************************************************************************************/
12321295

12331296
/********************************************************************************************/
@@ -1442,6 +1505,6 @@ std::string RawDataManager::describeEvent()
14421505
out << "## TF:Event " << mTimeFrameNo << ":" << mEventNo << ": "
14431506
// << hits->getsize() << " hits "
14441507
<< mTriggerRecord.getNumberOfDigits() << " digits and "
1445-
<< mTriggerRecord.getNumberOfTracklets() << " tracklets";
1508+
<< mTriggerRecord.getNumberOiTracklets() << " tracklets";
14461509
return out.str();
14471510
}

0 commit comments

Comments
 (0)