]>
Commit | Line | Data |
---|---|---|
c7585a2a | 1 | // $Id$ |
2 | ||
3 | //************************************************************************** | |
4 | //* This file is property of and copyright by the ALICE HLT Project * | |
5 | //* ALICE Experiment at CERN, All rights reserved. * | |
6 | //* * | |
7 | //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no> * | |
8 | //* for The ALICE HLT Project. * | |
9 | //* * | |
10 | //* Permission to use, copy, modify and distribute this software and its * | |
11 | //* documentation strictly for non-commercial purposes is hereby granted * | |
12 | //* without fee, provided that the above copyright notice appears in all * | |
13 | //* copies and that both the copyright notice and this permission notice * | |
14 | //* appear in the supporting documentation. The authors make no claims * | |
15 | //* about the suitability of this software for any purpose. It is * | |
16 | //* provided "as is" without express or implied warranty. * | |
17 | //************************************************************************** | |
18 | ||
19 | /// @file AliHLTTPCTrackGeometry.cxx | |
20 | /// @author Matthias Richter | |
21 | /// @date 2011-05-20 | |
22 | /// @brief Desciption of a track by a sequence of track points | |
23 | /// | |
24 | ||
25 | #include "AliHLTTPCTrackGeometry.h" | |
26 | #include "AliHLTTPCTransform.h" | |
27 | #include "AliHLTTPCSpacePointData.h" | |
eb50ad6b | 28 | #include "AliHLTTPCClusterDataFormat.h" |
29 | #include "AliHLTTPCSpacePointContainer.h" | |
b97434b7 | 30 | #include "AliHLTTPCHWCFSpacePointContainer.h" |
eb50ad6b | 31 | #include "AliHLTTPCDefinitions.h" |
32 | #include "AliHLTComponent.h" | |
c7585a2a | 33 | #include "AliHLTGlobalBarrelTrack.h" |
b97434b7 | 34 | #include "AliHLTDataDeflater.h" |
9b72add5 | 35 | #include "AliHLTErrorGuard.h" |
c7585a2a | 36 | #include "TMath.h" |
bbe2c314 | 37 | #include "TH2F.h" |
eb50ad6b | 38 | #include <memory> |
4a03ac27 | 39 | #include <sstream> |
40 | using namespace std; | |
c7585a2a | 41 | |
42 | /** ROOT macro for the implementation of ROOT specific class methods */ | |
43 | ClassImp(AliHLTTPCTrackGeometry) | |
44 | ||
45 | AliHLTTPCTrackGeometry::AliHLTTPCTrackGeometry() | |
46 | : AliHLTTrackGeometry() | |
bbe2c314 | 47 | , fRawTrackPoints() |
9b72add5 | 48 | , fDriftTimeFactorA(0.) |
49 | , fDriftTimeOffsetA(0.) | |
50 | , fDriftTimeFactorC(0.) | |
51 | , fDriftTimeOffsetC(0.) | |
c7585a2a | 52 | { |
53 | /// standard constructor | |
54 | } | |
55 | ||
56 | AliHLTTPCTrackGeometry::AliHLTTPCTrackGeometry(const AliHLTTPCTrackGeometry& src) | |
57 | : AliHLTTrackGeometry(src) | |
bbe2c314 | 58 | , fRawTrackPoints(src.fRawTrackPoints) |
9b72add5 | 59 | , fDriftTimeFactorA(0.) |
60 | , fDriftTimeOffsetA(0.) | |
61 | , fDriftTimeFactorC(0.) | |
62 | , fDriftTimeOffsetC(0.) | |
c7585a2a | 63 | { |
64 | /// copy constructor | |
65 | } | |
66 | ||
67 | AliHLTTPCTrackGeometry& AliHLTTPCTrackGeometry::operator=(const AliHLTTPCTrackGeometry& src) | |
68 | { | |
69 | /// assignment operator | |
70 | AliHLTTrackGeometry::operator=(src); | |
bbe2c314 | 71 | fRawTrackPoints.assign(src.fRawTrackPoints.begin(), src.fRawTrackPoints.end()); |
c7585a2a | 72 | return *this; |
73 | } | |
74 | ||
75 | AliHLTTPCTrackGeometry::~AliHLTTPCTrackGeometry() | |
76 | { | |
77 | /// destructor | |
78 | } | |
79 | ||
80 | float AliHLTTPCTrackGeometry::GetPlaneAlpha(AliHLTUInt32_t planeId) const | |
81 | { | |
82 | /// alpha of the plane | |
83 | UInt_t slice=AliHLTTPCSpacePointData::GetSlice(planeId); | |
84 | float alpha=( slice + 0.5 ) * TMath::Pi() / 9.0; | |
85 | if (alpha>TMath::TwoPi()) alpha-=TMath::TwoPi(); | |
86 | return alpha; | |
87 | } | |
88 | ||
89 | float AliHLTTPCTrackGeometry::GetPlaneR(AliHLTUInt32_t planeId) const | |
90 | { | |
91 | /// radial distance from global {0,0,0} | |
92 | UInt_t partition=AliHLTTPCSpacePointData::GetPatch(planeId); | |
93 | UInt_t number=AliHLTTPCSpacePointData::GetNumber(planeId); | |
94 | Int_t row=AliHLTTPCTransform::GetFirstRow(partition)+number; | |
95 | return AliHLTTPCTransform::Row2X(row); | |
96 | } | |
97 | ||
98 | float AliHLTTPCTrackGeometry::GetPlaneTheta(AliHLTUInt32_t /*planeId*/) const | |
99 | { | |
100 | /// theta of the plane | |
101 | return 0.0; | |
102 | } | |
103 | ||
104 | bool AliHLTTPCTrackGeometry::CheckBounds(AliHLTUInt32_t planeId, float u, float /*v*/) const | |
105 | { | |
106 | /// check bounds in u and v coordinate | |
107 | float r=GetPlaneR(planeId); | |
108 | if (r<AliHLTTPCTransform::GetFirstRow(0)) return false; | |
109 | ||
110 | // TODO: check if the pad width needs to be considered here | |
111 | return TMath::Abs(TMath::ASin(u/r))<=TMath::Pi()/18; | |
112 | } | |
113 | ||
114 | int AliHLTTPCTrackGeometry::CalculateTrackPoints(const AliHLTExternalTrackParam& track) | |
115 | { | |
116 | /// calculate the track points, expects the global magnetic field to be initialized | |
117 | AliHLTGlobalBarrelTrack bt(track); | |
118 | return CalculateTrackPoints(bt); | |
119 | } | |
120 | ||
121 | int AliHLTTPCTrackGeometry::CalculateTrackPoints(AliHLTGlobalBarrelTrack& track) | |
122 | { | |
123 | /// calculate the track points, expects the global magnetic field to be initialized | |
d5cf9283 | 124 | /// depending on the x coordinate of the first track point the corresponding padrow |
125 | /// is searched and the points are calculated outwards. Eventually the points are | |
126 | /// calculated inwards as a second step. | |
c7585a2a | 127 | int iResult=0; |
128 | int firstpadrow=0; | |
eb50ad6b | 129 | for (; |
130 | firstpadrow<AliHLTTPCTransform::GetNRows() && | |
131 | AliHLTTPCTransform::Row2X(firstpadrow)+AliHLTTPCTransform::GetPadLength(firstpadrow)<track.GetX(); | |
132 | firstpadrow++); | |
133 | if (firstpadrow>=AliHLTTPCTransform::GetNRows()) return 0; | |
d5cf9283 | 134 | // first calculated outwards |
c7585a2a | 135 | iResult=CalculateTrackPoints(track, firstpadrow, 1); |
136 | if (iResult>=0 && firstpadrow>0) | |
d5cf9283 | 137 | // now calculate inwards |
c7585a2a | 138 | iResult=CalculateTrackPoints(track, firstpadrow-1, -1); |
139 | return iResult; | |
140 | } | |
141 | ||
142 | int AliHLTTPCTrackGeometry::CalculateTrackPoints(AliHLTGlobalBarrelTrack& track, int firstpadrow, int step) | |
143 | { | |
144 | /// calculate the track points, expects the global magnetic field to be initialized | |
145 | float offsetAlpha=0.0; | |
146 | for (int padrow=firstpadrow; padrow>=0 && padrow<AliHLTTPCTransform::GetNRows(); padrow+=step) { | |
147 | float x=AliHLTTPCTransform::Row2X(padrow); | |
148 | float y=0.0; | |
149 | float z=0.0; | |
150 | ||
151 | int maxshift=9; | |
152 | int shift=0; | |
153 | int result=0; | |
154 | do { | |
155 | // start calculation of crossing points with padrow planes in the slice of the first point | |
156 | // plane alpha corresponds to alpha of the track, switch to neighboring slice if the result | |
157 | // is out of bounds | |
158 | if ((result=track.CalculateCrossingPoint(x, track.GetAlpha()-offsetAlpha, y, z))<1) break; | |
159 | float pointAlpha=TMath::ATan(y/x); | |
160 | if (TMath::Abs(pointAlpha)>TMath::Pi()/18) { | |
161 | offsetAlpha+=(pointAlpha>0?-1:1)*TMath::Pi()/9; | |
162 | result=0; | |
163 | } | |
164 | } while (result==0 && shift++<maxshift); | |
165 | if (result<1) continue; | |
166 | float planealpha=track.GetAlpha()-offsetAlpha; | |
167 | if (planealpha<0) planealpha+=TMath::TwoPi(); | |
6e65c8a6 | 168 | if (planealpha>TMath::TwoPi()) planealpha-=TMath::TwoPi(); |
102d9fef | 169 | int slice=int(9*planealpha/TMath::Pi()); |
b97434b7 | 170 | if (z<0) slice+=18; |
6e65c8a6 | 171 | if (slice>=36) { |
172 | HLTError("invalid slice %d calculated from alpha %f", slice, track.GetAlpha()); | |
173 | } | |
c7585a2a | 174 | int partition=AliHLTTPCTransform::GetPatch(padrow); |
175 | int row=padrow-AliHLTTPCTransform::GetFirstRow(partition); | |
176 | UInt_t id=AliHLTTPCSpacePointData::GetID(slice, partition, row); | |
177 | if (TMath::Abs(planealpha-GetPlaneAlpha(id))>0.0001) { | |
6e65c8a6 | 178 | HLTError("alpha missmatch for plane %08x (slice %d): alpha from id %f (%.0f deg), expected %f (%.0f deg)", id, slice, GetPlaneAlpha(id), 180*GetPlaneAlpha(id)/TMath::Pi(), planealpha, 180*planealpha/TMath::Pi()); |
c7585a2a | 179 | } |
bbe2c314 | 180 | if (AddTrackPoint(AliHLTTrackPoint(id, y, z), AliHLTTPCSpacePointData::GetID(slice, partition, 0))>=0) { |
181 | Float_t rpt[3]={0.,y,z}; // row pad time | |
182 | AliHLTTPCTransform::LocHLT2Raw(rpt, slice, padrow); | |
9b72add5 | 183 | float m=fDriftTimeFactorA; |
184 | float n=fDriftTimeOffsetA; | |
185 | if (slice>=18) { | |
186 | m=fDriftTimeFactorC; | |
187 | n=fDriftTimeOffsetC; | |
188 | } | |
189 | if (TMath::Abs(m)>0.) { | |
190 | rpt[2]=(z-n)/m; | |
4a03ac27 | 191 | if (step>0) { |
192 | // ascending padrows added at end | |
193 | fRawTrackPoints.push_back(AliHLTTrackPoint(id, rpt[1], rpt[2])); | |
194 | } else { | |
195 | // descending padrows added at begin | |
196 | fRawTrackPoints.insert(fRawTrackPoints.begin(), AliHLTTrackPoint(id, rpt[1], rpt[2])); | |
197 | } | |
198 | // FIXME: implement a Print function for the raw track points | |
199 | // stringstream sout; | |
200 | // sout << " slice " << setfill(' ') << setw(3) << fixed << right << slice | |
201 | // << " row " << setfill(' ') << setw(3) << fixed << right << padrow | |
202 | // << " id 0x" << hex << setw(8) << id << dec | |
6e65c8a6 | 203 | // << " y " << setfill(' ') << setw(5) << fixed << right << setprecision (2) << y |
204 | // << " z " << setfill(' ') << setw(5) << fixed << right << setprecision (2) << z | |
205 | // << " pad " << setfill(' ') << setw(5) << fixed << right << setprecision (2) << rpt[1] | |
206 | // << " time " << setfill(' ') << setw(5) << fixed << right << setprecision (2) << rpt[2] | |
207 | // << endl; | |
4a03ac27 | 208 | // cout << sout.str(); |
9b72add5 | 209 | } else { |
210 | ALIHLTERRORGUARD(1, "drift time correction not initialized, can not add track points in raw coordinates"); | |
211 | } | |
bbe2c314 | 212 | } |
c7585a2a | 213 | } |
214 | return 0; | |
215 | } | |
216 | ||
217 | int AliHLTTPCTrackGeometry::FindMatchingTrackPoint(AliHLTUInt32_t spacepointId, float spacepoint[3], AliHLTUInt32_t& planeId) | |
218 | { | |
219 | /// find the track point which can be associated to a spacepoint with coordinates and id | |
220 | UInt_t slice=AliHLTTPCSpacePointData::GetSlice(spacepointId); | |
221 | UInt_t partition=AliHLTTPCSpacePointData::GetPatch(spacepointId); | |
222 | int row=AliHLTTPCTransform::GetPadRow(spacepoint[0]); | |
9003c201 | 223 | bool bSpecialRow=row==30 || row==90 || row==139; |
c7585a2a | 224 | if (row<AliHLTTPCTransform::GetFirstRow(partition) || row>AliHLTTPCTransform::GetLastRow(partition)) { |
225 | HLTError("row number %d calculated from x value %f is outside slice %d partition %d", row, spacepoint[0], slice, partition); | |
226 | return -EINVAL; | |
227 | } | |
9003c201 | 228 | |
229 | // find the crossing point of the track with the padrow plane where | |
230 | // the spacepoint is | |
231 | // 1) calculate plane id from slice, partition and row (within partition) | |
c7585a2a | 232 | row-=AliHLTTPCTransform::GetFirstRow(partition); |
233 | UInt_t id=AliHLTTPCSpacePointData::GetID(slice, partition, row); | |
234 | const AliHLTTrackPoint* point=GetTrackPoint(id); | |
9003c201 | 235 | // track might be outside the partition and cross the central membrane |
236 | // search in the other half of the TPC | |
c7585a2a | 237 | if (!point && slice<18) { |
9003c201 | 238 | // search in the neighboring partition on the C side |
c7585a2a | 239 | id=AliHLTTPCSpacePointData::GetID(slice+18, partition, row); |
240 | point=GetTrackPoint(id); | |
241 | } else if (!point && slice>=18) { | |
9003c201 | 242 | // search in the neighboring partition on the A side |
c7585a2a | 243 | id=AliHLTTPCSpacePointData::GetID(slice-18, partition, row); |
244 | point=GetTrackPoint(id); | |
245 | } | |
9003c201 | 246 | |
247 | // search in the neighboring partition, this takes account for rows | |
248 | // 30, 90, and 139 which are partly in one and the other partition | |
249 | if (!point && bSpecialRow) { | |
250 | row+=AliHLTTPCTransform::GetFirstRow(partition); | |
251 | row-=AliHLTTPCTransform::GetFirstRow(partition-1); | |
252 | id=AliHLTTPCSpacePointData::GetID(slice, partition-1, row); | |
253 | point=GetTrackPoint(id); | |
254 | if (!point && slice<18) { | |
255 | // search in the neighboring partition on the C side | |
256 | id=AliHLTTPCSpacePointData::GetID(slice+18, partition-1, row); | |
257 | point=GetTrackPoint(id); | |
258 | } else if (!point && slice>=18) { | |
259 | // search in the neighboring partition on the A side | |
260 | id=AliHLTTPCSpacePointData::GetID(slice-18, partition-1, row); | |
261 | point=GetTrackPoint(id); | |
262 | } | |
263 | } | |
264 | ||
c7585a2a | 265 | if (point) { |
266 | planeId=id; | |
9003c201 | 267 | if (point->HaveAssociatedSpacePoint()) { |
268 | if (GetVerbosity()>2) HLTInfo("descarding spacepoint 0x%08x z=%f y=%f z=%f: track point 0x%08x already occupied", spacepoint[0], spacepoint[1], spacepoint[2], planeId); | |
269 | return 0; // already occupied | |
270 | } | |
271 | float maxdy=2.; | |
272 | float maxdz=2.; | |
273 | if (TMath::Abs(point->GetU()-spacepoint[1])>maxdy) { | |
274 | if (GetVerbosity()>0) HLTInfo("descarding spacepoint 0x%08x y=%f z=%f: track point 0x%08x y %f outside tolerance %f", spacepoint[1], spacepoint[2], planeId, point->GetU(), maxdy); | |
275 | return -ENOENT; | |
276 | } | |
277 | if (TMath::Abs(point->GetV()-spacepoint[2])>maxdz) { | |
278 | if (GetVerbosity()>0) HLTInfo("descarding spacepoint 0x%08x y=%f z=%f: track point 0x%08x z %f outside tolerance %f", spacepoint[1], spacepoint[2], planeId, point->GetV(), maxdz); | |
279 | return -ENOENT; | |
280 | } | |
c7585a2a | 281 | return 1; |
282 | } | |
283 | return -ENOENT; | |
284 | } | |
eb50ad6b | 285 | |
c2bea8b3 | 286 | |
287 | int AliHLTTPCTrackGeometry::RegisterTrackPoints(AliHLTTrackGrid* pGrid) const | |
288 | { | |
289 | /// register track points in the index grid, at this step the number | |
290 | /// of tracks in each cell is counted | |
291 | if (!pGrid) return -EINVAL; | |
292 | int iResult=0; | |
293 | for (vector<AliHLTTrackPoint>::const_iterator tp=TrackPoints().begin(); | |
294 | tp!=TrackPoints().end() && iResult>=0; tp++) { | |
295 | AliHLTUInt32_t id=tp->GetId(); | |
296 | iResult=pGrid->CountSpacePoint(AliHLTTPCSpacePointData::GetSlice(id), | |
297 | AliHLTTPCSpacePointData::GetPatch(id), | |
298 | AliHLTTPCSpacePointData::GetNumber(id)); | |
299 | } | |
300 | return iResult; | |
301 | } | |
302 | ||
303 | int AliHLTTPCTrackGeometry::FillTrackPoints(AliHLTTrackGrid* pGrid) const | |
304 | { | |
305 | /// fill track points to index grid | |
306 | if (!pGrid) return -EINVAL; | |
307 | int iResult=0; | |
308 | for (vector<AliHLTTrackPoint>::const_iterator tp=TrackPoints().begin(); | |
309 | tp!=TrackPoints().end() && iResult>=0; tp++) { | |
310 | AliHLTUInt32_t id=tp->GetId(); | |
311 | iResult=pGrid->AddSpacePoint(GetTrackId(), | |
312 | AliHLTTPCSpacePointData::GetSlice(id), | |
313 | AliHLTTPCSpacePointData::GetPatch(id), | |
314 | AliHLTTPCSpacePointData::GetNumber(id)); | |
315 | } | |
316 | return iResult; | |
317 | } | |
318 | ||
9003c201 | 319 | AliHLTSpacePointContainer* AliHLTTPCTrackGeometry::ConvertToSpacePoints(bool bAssociated) const |
eb50ad6b | 320 | { |
321 | /// create a collection of all points | |
322 | std::auto_ptr<AliHLTTPCSpacePointContainer> spacepoints(new AliHLTTPCSpacePointContainer); | |
323 | if (!spacepoints.get()) return NULL; | |
324 | ||
dd2011bf | 325 | const vector<AliHLTTrackPoint>& trackPoints=TrackPoints(); |
eb50ad6b | 326 | unsigned i=0; |
327 | while (i<trackPoints.size()) { | |
328 | // allocate buffer for all points, even though the buffer might not be filled | |
329 | // completely because of a partition change | |
330 | int nofPoints=trackPoints.size()-i; | |
331 | int blocksize=sizeof(AliHLTTPCClusterData)+nofPoints*sizeof(AliHLTTPCSpacePointData); | |
332 | AliHLTUInt8_t* pBuffer=spacepoints->Alloc(blocksize); | |
333 | if (!pBuffer) return NULL; | |
334 | AliHLTTPCClusterData* pClusterData=reinterpret_cast<AliHLTTPCClusterData*>(pBuffer); | |
335 | pClusterData->fSpacePointCnt=0; | |
336 | AliHLTTPCSpacePointData* pClusters=pClusterData->fSpacePoints; | |
337 | int currentSlice=-1; | |
338 | int currentPartition=-1; | |
339 | for (; i<trackPoints.size(); i++) { | |
9003c201 | 340 | if (bAssociated && !trackPoints[i].HaveAssociatedSpacePoint()) continue; |
eb50ad6b | 341 | AliHLTUInt32_t planeId=trackPoints[i].GetId(); |
342 | int slice=AliHLTTPCSpacePointData::GetSlice(planeId); | |
343 | int partition=AliHLTTPCSpacePointData::GetPatch(planeId); | |
344 | int number=AliHLTTPCSpacePointData::GetNumber(planeId); | |
345 | if ((currentSlice>=0 && currentSlice!=slice) || (currentPartition>=0 && currentPartition!=partition)) { | |
346 | // change of partition or slice, need to go to next block | |
347 | // 2011-07-26 currently all spacepoints go into one block, if separated | |
348 | // blocks per partition are needed one has to leave the inner loop here | |
349 | // and set the data block specification below | |
350 | // Caution: not tested, only the last block seems to make it through | |
351 | //break; | |
352 | } | |
353 | currentSlice=slice; | |
354 | currentPartition=partition; | |
355 | pClusters[pClusterData->fSpacePointCnt].fX=GetPlaneR(planeId); | |
356 | pClusters[pClusterData->fSpacePointCnt].fY=trackPoints[i].GetU(); | |
357 | pClusters[pClusterData->fSpacePointCnt].fZ=trackPoints[i].GetV(); | |
358 | pClusters[pClusterData->fSpacePointCnt].fID=planeId; | |
359 | pClusters[pClusterData->fSpacePointCnt].fPadRow=AliHLTTPCTransform::GetFirstRow(partition)+number; | |
360 | pClusters[pClusterData->fSpacePointCnt].fSigmaY2=0.; | |
361 | pClusters[pClusterData->fSpacePointCnt].fSigmaZ2=0.; | |
362 | pClusters[pClusterData->fSpacePointCnt].fCharge=0; | |
363 | pClusters[pClusterData->fSpacePointCnt].fQMax=0; | |
364 | pClusters[pClusterData->fSpacePointCnt].fUsed=0; | |
365 | pClusters[pClusterData->fSpacePointCnt].fTrackN=0; | |
366 | pClusterData->fSpacePointCnt++; | |
367 | } | |
368 | AliHLTComponentBlockData bd; | |
369 | AliHLTComponent::FillBlockData(bd); | |
370 | bd.fPtr=pBuffer; | |
9003c201 | 371 | bd.fSize=sizeof(AliHLTTPCClusterData)+pClusterData->fSpacePointCnt*sizeof(AliHLTTPCSpacePointData); |
eb50ad6b | 372 | AliHLTComponent::SetDataType(bd.fDataType, "CLUSTERS", "TPC "); |
bbe2c314 | 373 | bd.fSpecification=kAliHLTVoidDataSpec;//AliHLTTPCDefinitions::EncodeDataSpecification(currentSlice, currentSlice, currentPartition, currentPartition); |
eb50ad6b | 374 | spacepoints->AddInputBlock(&bd); |
375 | } | |
376 | ||
377 | return spacepoints.release(); | |
378 | } | |
bbe2c314 | 379 | |
c2bea8b3 | 380 | const AliHLTTrackGeometry::AliHLTTrackPoint* AliHLTTPCTrackGeometry::GetRawTrackPoint(AliHLTUInt32_t id) const |
381 | { | |
382 | /// get raw track point of id | |
383 | const AliHLTTrackGeometry::AliHLTTrackPoint* p=find(&fRawTrackPoints[0], &fRawTrackPoints[fRawTrackPoints.size()], id); | |
384 | if (p==&fRawTrackPoints[fRawTrackPoints.size()]) return 0; | |
385 | return p; | |
386 | } | |
387 | ||
388 | AliHLTTrackGeometry::AliHLTTrackPoint* AliHLTTPCTrackGeometry::GetRawTrackPoint(AliHLTUInt32_t id) | |
389 | { | |
390 | /// get raw track point of id | |
391 | AliHLTTrackGeometry::AliHLTTrackPoint* p=find(&fRawTrackPoints[0], &fRawTrackPoints[fRawTrackPoints.size()], id); | |
392 | if (p==&fRawTrackPoints[fRawTrackPoints.size()]) return 0; | |
393 | return p; | |
394 | } | |
395 | ||
bbe2c314 | 396 | int AliHLTTPCTrackGeometry::FillRawResidual(int coordinate, TH2* histo, AliHLTSpacePointContainer* points) const |
397 | { | |
398 | // fill residual histogram | |
399 | if (!histo || !points) return -EINVAL; | |
400 | const vector<AliHLTTrackPoint>& trackPoints=TrackPoints(); | |
401 | for (vector<AliHLTTrackPoint>::const_iterator trackpoint=trackPoints.begin(); | |
402 | trackpoint!=trackPoints.end(); trackpoint++) { | |
6c6615c1 | 403 | if (!trackpoint->HaveAssociatedSpacePoint()) continue; |
404 | for (vector<AliHLTTrackSpacepoint>::const_iterator sp=(trackpoint->GetSpacepoints()).begin(); | |
405 | sp!=(trackpoint->GetSpacepoints()).end(); sp++) { | |
406 | AliHLTUInt32_t spacepointId=sp->fId; | |
bbe2c314 | 407 | vector<AliHLTTrackPoint>::const_iterator rawpoint=find(fRawTrackPoints.begin(), fRawTrackPoints.end(), trackpoint->GetId()); |
408 | if (rawpoint==fRawTrackPoints.end()) { | |
409 | HLTError("can not find track raw coordinates of track point 0x%08x", trackpoint->GetId()); | |
410 | continue; | |
411 | } | |
412 | if (!points->Check(spacepointId)) { | |
413 | //HLTError("can not find associated space point 0x%08x of track point 0x%08x", spacepointId, trackpoint->GetId()); | |
414 | continue; | |
415 | } | |
416 | float value=0.; | |
417 | if (coordinate==0) { | |
418 | value=rawpoint->GetU()-points->GetY(spacepointId); | |
419 | histo->Fill(GetPlaneR(trackpoint->GetId()), value); | |
420 | } else { | |
421 | value=rawpoint->GetV()-points->GetZ(spacepointId); | |
422 | //histo->Fill(GetPlaneR(trackpoint->GetId()), value); | |
423 | histo->Fill(rawpoint->GetV(), value); | |
424 | } | |
6c6615c1 | 425 | } |
bbe2c314 | 426 | } |
427 | return 0; | |
428 | } | |
b97434b7 | 429 | |
430 | int AliHLTTPCTrackGeometry::Write(const AliHLTGlobalBarrelTrack& track, | |
431 | AliHLTSpacePointContainer* pSpacePoints, | |
432 | AliHLTDataDeflater* pDeflater, | |
433 | AliHLTUInt8_t* outputPtr, | |
434 | AliHLTUInt32_t size, | |
9b72add5 | 435 | vector<AliHLTUInt32_t>* writtenClusterIds, |
b97434b7 | 436 | const char* option) const |
437 | { | |
438 | // write track block to buffer | |
439 | if (size<=sizeof(AliHLTTPCTrackBlock)) return -ENOSPC; | |
440 | AliHLTTPCTrackBlock* pTrackBlock=reinterpret_cast<AliHLTTPCTrackBlock*>(outputPtr); | |
441 | pTrackBlock->fSize=sizeof(AliHLTTPCTrackBlock); // size of cluster block added later | |
6e65c8a6 | 442 | float alpha=track.GetAlpha(); |
443 | while (alpha<0.) alpha+=TMath::TwoPi(); | |
444 | while (alpha>TMath::TwoPi()) alpha-=TMath::TwoPi(); | |
445 | pTrackBlock->fSlice=AliHLTUInt8_t(9*alpha/TMath::Pi()); | |
446 | if (pTrackBlock->fSlice>=36) { | |
447 | HLTError("invalid slice %d calculated from alpha %f", pTrackBlock->fSlice, track.GetAlpha()); | |
448 | } | |
b97434b7 | 449 | pTrackBlock->fReserved=0; |
450 | pTrackBlock->fX = track.GetX(); | |
451 | pTrackBlock->fY = track.GetY(); | |
452 | pTrackBlock->fZ = track.GetZ(); | |
453 | pTrackBlock->fSinPsi = track.GetSnp(); | |
454 | pTrackBlock->fTgl = track.GetTgl(); | |
455 | pTrackBlock->fq1Pt = track.GetSigned1Pt(); | |
456 | ||
457 | pDeflater->Clear(); | |
458 | pDeflater->InitBitDataOutput(reinterpret_cast<AliHLTUInt8_t*>(outputPtr+sizeof(AliHLTTPCTrackBlock)), size-sizeof(AliHLTTPCTrackBlock)); | |
9b72add5 | 459 | int result=WriteAssociatedClusters(pSpacePoints, pDeflater, writtenClusterIds, option); |
b97434b7 | 460 | if (result<0) return result; |
461 | pTrackBlock->fSize+=result; | |
462 | return pTrackBlock->fSize; | |
463 | } | |
464 | ||
465 | int AliHLTTPCTrackGeometry::WriteAssociatedClusters(AliHLTSpacePointContainer* pSpacePoints, | |
466 | AliHLTDataDeflater* pDeflater, | |
9b72add5 | 467 | vector<AliHLTUInt32_t>* writtenClusterIds, |
b97434b7 | 468 | const char* /*option*/) const |
469 | { | |
470 | // write associated clusters to buffer via deflater | |
471 | if (!pDeflater || !pSpacePoints) return -EINVAL; | |
472 | AliHLTTPCHWCFSpacePointContainer* pTPCRawSpacePoints=dynamic_cast<AliHLTTPCHWCFSpacePointContainer*>(pSpacePoints); | |
473 | if (!pTPCRawSpacePoints) return -EINVAL; | |
474 | bool bReverse=true; | |
475 | bool bWriteSuccess=true; | |
476 | int writtenClusters=0; | |
477 | // filling of track points starts from first point on track outwards, and | |
478 | // then from that point inwards. That's why the lower padrows might be in | |
479 | // reverse order at the end of the track point array. If the last element | |
480 | // is bigger than the first element, only trackpoints in ascending order | |
481 | // are in the array | |
482 | vector<AliHLTTrackPoint>::const_iterator clrow=fRawTrackPoints.end(); | |
483 | if (clrow!=fRawTrackPoints.begin()) { | |
484 | clrow--; | |
485 | AliHLTUInt32_t partition=AliHLTTPCSpacePointData::GetPatch(clrow->GetId()); | |
486 | AliHLTUInt32_t partitionrow=AliHLTTPCSpacePointData::GetNumber(clrow->GetId()); | |
487 | partitionrow+=AliHLTTPCTransform::GetFirstRow(partition); | |
488 | AliHLTUInt32_t firstpartition=AliHLTTPCSpacePointData::GetPatch(fRawTrackPoints.begin()->GetId()); | |
489 | AliHLTUInt32_t firstpartitionrow=AliHLTTPCSpacePointData::GetNumber(fRawTrackPoints.begin()->GetId()); | |
490 | firstpartitionrow+=AliHLTTPCTransform::GetFirstRow(firstpartition); | |
491 | if (partitionrow>=firstpartitionrow) { | |
492 | bReverse=false; | |
493 | clrow=fRawTrackPoints.begin(); | |
494 | } | |
495 | } | |
6e65c8a6 | 496 | const AliHLTUInt32_t clusterCountBitLength=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[AliHLTTPCDefinitions::kClusterCount].fBitLength; |
b97434b7 | 497 | unsigned long dataPosition=pDeflater->GetCurrentByteOutputPosition(); |
498 | for (unsigned row=0; row<159 && bWriteSuccess; row++) { | |
499 | if (clrow!=fRawTrackPoints.end()) { | |
500 | AliHLTUInt32_t thisPartition=AliHLTTPCSpacePointData::GetPatch(clrow->GetId()); | |
501 | AliHLTUInt32_t thisTrackRow=AliHLTTPCSpacePointData::GetNumber(clrow->GetId()); | |
502 | thisTrackRow+=AliHLTTPCTransform::GetFirstRow(thisPartition); | |
503 | if (thisTrackRow==row) { | |
504 | // write clusters | |
505 | const vector<AliHLTTrackSpacepoint>& clusters=clrow->GetSpacepoints(); | |
506 | AliHLTUInt32_t haveClusters=clusters.size()>0; | |
507 | // 1 bit for clusters on that padrow | |
508 | bWriteSuccess=bWriteSuccess && pDeflater->OutputBit(haveClusters); | |
509 | if (haveClusters) { | |
6e65c8a6 | 510 | bWriteSuccess=bWriteSuccess && pDeflater->OutputBits(clusters.size(), clusterCountBitLength); |
b97434b7 | 511 | for (vector<AliHLTTrackSpacepoint>::const_iterator clid=clusters.begin(); |
512 | clid!=clusters.end() && bWriteSuccess; clid++) { | |
513 | if (!pSpacePoints->Check(clid->fId)) { | |
514 | HLTError("can not find spacepoint 0x%08x", clid->fId); | |
515 | continue; | |
516 | } | |
9b72add5 | 517 | if (writtenClusterIds) { |
518 | writtenClusterIds->push_back(clid->fId); | |
519 | } | |
b97434b7 | 520 | |
6e65c8a6 | 521 | // FIXME: there is a bug in the calculation of the residuals stored with the |
522 | // assiciated space point, calculate again, but needs to be fixed | |
523 | float deltapad =pSpacePoints->GetY(clid->fId)-clrow->GetU();//clid->fdU; | |
524 | float deltatime =pSpacePoints->GetZ(clid->fId)-clrow->GetV();//clid->fdV; | |
4a03ac27 | 525 | if (TMath::Abs(deltapad)>=AliHLTTPCDefinitions::fgkMaxClusterDeltaPad || |
526 | TMath::Abs(deltatime)>=AliHLTTPCDefinitions::fgkMaxClusterDeltaTime) { | |
527 | AliHLTUInt8_t slice = AliHLTTPCSpacePointData::GetSlice(clid->fId); | |
528 | AliHLTUInt8_t partition = AliHLTTPCSpacePointData::GetPatch(clid->fId); | |
529 | HLTFatal("cluster 0x%08x slice %d partition %d: residual out of range - pad %f (max %d), time %f (max %d)", clid->fId, slice, partition, deltapad, AliHLTTPCDefinitions::fgkMaxClusterDeltaPad, deltatime, AliHLTTPCDefinitions::fgkMaxClusterDeltaTime); | |
530 | } | |
b97434b7 | 531 | float sigmaY2=pSpacePoints->GetYWidth(clid->fId); |
532 | float sigmaZ2=pSpacePoints->GetZWidth(clid->fId); | |
533 | AliHLTUInt64_t charge=(AliHLTUInt64_t)pSpacePoints->GetCharge(clid->fId); | |
534 | AliHLTUInt64_t qmax=(AliHLTUInt64_t)pTPCRawSpacePoints->GetQMax(clid->fId); | |
6e65c8a6 | 535 | // cout << " row " << setfill(' ') << setw(3) << fixed << right << row |
536 | // << " pad " << setfill(' ') << setw(7) << fixed << right << setprecision (4) << pSpacePoints->GetY(clid->fId) | |
537 | // << " dpad " << setfill(' ') << setw(7) << fixed << right << setprecision (4) << deltapad | |
538 | // << " time " << setfill(' ') << setw(7) << fixed << right << setprecision (4) << pSpacePoints->GetZ(clid->fId) | |
539 | // << " dtime " << setfill(' ') << setw(7) << fixed << right << setprecision (4) << deltatime | |
540 | // << " charge " << setfill(' ') << setw(5) << fixed << right << charge | |
541 | // << " qmax " << setfill(' ') << setw(4) << fixed << right << qmax | |
542 | // << endl; | |
b97434b7 | 543 | |
8a3426fd | 544 | // time and pad coordinates are scaled and transformed to integer values for |
545 | // both cluster and track point before calculating the residual. this makes | |
546 | // the compression lossless with respect to the format without track model | |
547 | // compression | |
b97434b7 | 548 | AliHLTUInt64_t deltapad64=0; |
549 | AliHLTUInt32_t signDeltaPad=0; | |
550 | if (!isnan(deltapad)) { | |
8a3426fd | 551 | double clusterpad=pSpacePoints->GetY(clid->fId); |
552 | double trackpad=clrow->GetU(); | |
553 | if (clusterpad<0.) { | |
554 | HLTError("cluster 0x%08x has negative pad position", clid->fId, clusterpad); | |
555 | } | |
556 | clusterpad*=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[AliHLTTPCDefinitions::kResidualPad].fScale; | |
557 | trackpad*=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[AliHLTTPCDefinitions::kResidualPad].fScale; | |
558 | AliHLTUInt64_t clusterpad64=(AliHLTUInt64_t)round(clusterpad); | |
559 | AliHLTUInt64_t trackpad64=0; | |
560 | if (trackpad>0.) trackpad64=(AliHLTUInt64_t)round(trackpad); | |
561 | if (clusterpad64<trackpad64) { | |
562 | deltapad64=trackpad64-clusterpad64; | |
563 | signDeltaPad=1; | |
564 | } else { | |
565 | deltapad64=clusterpad64-trackpad64; | |
566 | signDeltaPad=0; | |
567 | } | |
b97434b7 | 568 | } |
569 | AliHLTUInt64_t deltatime64=0; | |
570 | AliHLTUInt32_t signDeltaTime=0; | |
571 | if (!isnan(deltatime)) { | |
8a3426fd | 572 | double clustertime=pSpacePoints->GetZ(clid->fId); |
573 | double tracktime=clrow->GetV(); | |
574 | if (clustertime<0.) { | |
575 | HLTError("cluster 0x%08x has negative time position", clid->fId, clustertime); | |
576 | } | |
577 | clustertime*=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[AliHLTTPCDefinitions::kResidualTime].fScale; | |
578 | tracktime*=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[AliHLTTPCDefinitions::kResidualTime].fScale; | |
579 | AliHLTUInt64_t clustertime64=(AliHLTUInt64_t)round(clustertime); | |
580 | AliHLTUInt64_t tracktime64=0; | |
581 | if (tracktime>0.) tracktime64=(AliHLTUInt64_t)round(tracktime); | |
582 | if (clustertime64<tracktime64) { | |
583 | deltatime64=tracktime64-clustertime64; | |
584 | signDeltaTime=1; | |
585 | } else { | |
586 | deltatime64=clustertime64-tracktime64; | |
587 | signDeltaTime=0; | |
588 | } | |
b97434b7 | 589 | } |
590 | AliHLTUInt64_t sigmaY264=0; | |
591 | if (!isnan(sigmaY2)) sigmaY264=(AliHLTUInt64_t)round(sigmaY2*AliHLTTPCDefinitions::fgkClusterParameterDefinitions[AliHLTTPCDefinitions::kSigmaY2].fScale); | |
592 | AliHLTUInt64_t sigmaZ264=0; | |
593 | if (!isnan(sigmaZ2)) sigmaZ264=(AliHLTUInt64_t)round(sigmaZ2*AliHLTTPCDefinitions::fgkClusterParameterDefinitions[AliHLTTPCDefinitions::kSigmaZ2].fScale); | |
b97434b7 | 594 | bWriteSuccess=bWriteSuccess && pDeflater->OutputParameterBits(AliHLTTPCDefinitions::kResidualPad , deltapad64); |
6e65c8a6 | 595 | bWriteSuccess=bWriteSuccess && pDeflater->OutputBit(signDeltaPad); |
b97434b7 | 596 | bWriteSuccess=bWriteSuccess && pDeflater->OutputParameterBits(AliHLTTPCDefinitions::kResidualTime , deltatime64); |
6e65c8a6 | 597 | bWriteSuccess=bWriteSuccess && pDeflater->OutputBit(signDeltaTime); |
b97434b7 | 598 | bWriteSuccess=bWriteSuccess && pDeflater->OutputParameterBits(AliHLTTPCDefinitions::kSigmaY2 , sigmaY264); |
599 | bWriteSuccess=bWriteSuccess && pDeflater->OutputParameterBits(AliHLTTPCDefinitions::kSigmaZ2 , sigmaZ264); | |
600 | bWriteSuccess=bWriteSuccess && pDeflater->OutputParameterBits(AliHLTTPCDefinitions::kCharge , charge); | |
601 | bWriteSuccess=bWriteSuccess && pDeflater->OutputParameterBits(AliHLTTPCDefinitions::kQMax , qmax); | |
602 | if (bWriteSuccess) writtenClusters++; | |
603 | } | |
604 | } | |
605 | ||
606 | // set to next trackpoint | |
607 | if (bReverse) { | |
608 | if (clrow!=fRawTrackPoints.begin()) { | |
609 | AliHLTUInt32_t nextPartition=AliHLTTPCSpacePointData::GetPatch((clrow-1)->GetId()); | |
610 | AliHLTUInt32_t nextTrackRow=AliHLTTPCSpacePointData::GetNumber((clrow-1)->GetId()); | |
611 | nextTrackRow+=AliHLTTPCTransform::GetFirstRow(nextPartition); | |
612 | if (thisTrackRow+1==nextTrackRow) { | |
613 | clrow--; | |
614 | } else { | |
615 | // switch direction start from beginning | |
616 | clrow=fRawTrackPoints.begin(); | |
617 | bReverse=false; | |
618 | } | |
619 | } else { | |
620 | // all trackpoints processed | |
621 | clrow=fRawTrackPoints.end(); | |
622 | } | |
623 | } else { | |
624 | clrow++; | |
625 | } | |
626 | continue; | |
627 | } else { | |
628 | // sequence not ordered, search | |
629 | // this has been fixed and the search is no longer necessary | |
630 | // for (clrow=fRawTrackPoints.begin(); clrow!=fRawTrackPoints.end(); clrow++) { | |
631 | // if ((AliHLTTPCSpacePointData::GetNumber(clrow->GetId())+AliHLTTPCTransform::GetFirstRow(AliHLTTPCSpacePointData::GetPatch(clrow->GetId())))==row) break; | |
632 | // } | |
633 | // if (clrow==fRawTrackPoints.end()) { | |
634 | // clrow=fRawTrackPoints.begin(); | |
635 | // HLTWarning("no trackpoint on row %d, current point %d", row, thisTrackRow); | |
636 | // } | |
637 | } | |
638 | } | |
639 | // no cluster on that padrow | |
640 | AliHLTUInt32_t haveClusters=0; | |
641 | bWriteSuccess=bWriteSuccess && pDeflater->OutputBit(haveClusters); | |
642 | } | |
643 | ||
644 | if (!bWriteSuccess) return -ENOSPC; | |
645 | ||
646 | int allClusters=0; | |
647 | for (clrow=fRawTrackPoints.begin(); clrow!=fRawTrackPoints.end(); clrow++) { | |
648 | allClusters+=clrow->GetSpacepoints().size(); | |
649 | } | |
650 | if (allClusters!=writtenClusters) { | |
651 | HLTError("track %d mismatch in written clusters: %d but expected %d", GetTrackId(), writtenClusters, allClusters); | |
652 | } | |
653 | ||
654 | pDeflater->Pad8Bits(); | |
655 | return pDeflater->GetCurrentByteOutputPosition()-dataPosition; | |
656 | } | |
657 | ||
6e65c8a6 | 658 | int AliHLTTPCTrackGeometry::Read(const AliHLTUInt8_t* buffer, |
659 | AliHLTUInt32_t size, | |
660 | float bz, | |
661 | AliHLTUInt32_t& clusterBlockSize, | |
662 | const char* /*option*/) | |
663 | { | |
664 | // read track block from buffer | |
665 | int iResult=0; | |
666 | if (!buffer) return -EINVAL; | |
667 | if (size<sizeof(AliHLTTPCTrackBlock)) { | |
668 | HLTError("buffer does not contain valid data of track model clusters"); | |
669 | return -ENODATA; | |
670 | } | |
671 | const AliHLTTPCTrackBlock* pTrackBlock=reinterpret_cast<const AliHLTTPCTrackBlock*>(buffer); | |
672 | if (pTrackBlock->fSize>size) { | |
673 | HLTError("inconsistent track data block of size %d exceeds available buffer of size %d", pTrackBlock->fSize, size); | |
674 | return -ENODATA; | |
675 | } | |
676 | if (pTrackBlock->fSize<sizeof(AliHLTTPCTrackBlock)) { | |
677 | HLTError("inconsistent size of track data block specified in the header: %d", pTrackBlock->fSize); | |
678 | return -ENODATA; | |
679 | } | |
680 | AliHLTExternalTrackParam param; | |
681 | memset(¶m, 0, sizeof(param)); | |
682 | param.fAlpha =( pTrackBlock->fSlice + 0.5 ) * TMath::Pi() / 9.0; | |
683 | if (param.fAlpha>TMath::TwoPi()) param.fAlpha-=TMath::TwoPi(); | |
684 | param.fX = pTrackBlock->fX; | |
685 | param.fY = pTrackBlock->fY; | |
686 | param.fZ = pTrackBlock->fZ; | |
687 | param.fSinPsi = pTrackBlock->fSinPsi; | |
688 | param.fTgl = pTrackBlock->fTgl; | |
689 | param.fq1Pt = pTrackBlock->fq1Pt; | |
690 | AliHLTGlobalBarrelTrack track(param); | |
691 | if ((iResult=track.CalculateHelixParams(bz))<0) { | |
692 | HLTError("failed to calculate helix params: %d", iResult); | |
693 | return iResult; | |
694 | } | |
695 | if ((iResult=CalculateTrackPoints(track))<0) { | |
696 | HLTError("failed to calculate track points: %d", iResult); | |
697 | return iResult; | |
698 | } | |
699 | clusterBlockSize=pTrackBlock->fSize-sizeof(AliHLTTPCTrackBlock); | |
700 | return sizeof(AliHLTTPCTrackBlock); | |
701 | } |