Line data Source code
1 : //**************************************************************************
2 : //* This file is property of and copyright by the ALICE HLT Project *
3 : //* ALICE Experiment at CERN, All rights reserved. *
4 : //* *
5 : //* Primary Authors: Sergey Gorbunov <sergey.gorbunov@cern.ch> *
6 : //* for The ALICE HLT Project. *
7 : //* *
8 : //* Permission to use, copy, modify and distribute this software and its *
9 : //* documentation strictly for non-commercial purposes is hereby granted *
10 : //* without fee, provided that the above copyright notice appears in all *
11 : //* copies and that both the copyright notice and this permission notice *
12 : //* appear in the supporting documentation. The authors make no claims *
13 : //* about the suitability of this software for any purpose. It is *
14 : //* provided "as is" without express or implied warranty. *
15 : //**************************************************************************
16 :
17 : /** @file AliHLTTPCFastTransform.cxx
18 : @author Sergey Gorbubnov
19 : @date
20 : @brief
21 : */
22 :
23 :
24 : #include "AliHLTTPCFastTransform.h"
25 : #include "AliTPCTransform.h"
26 : #include "AliTPCParam.h"
27 : #include "AliTPCRecoParam.h"
28 : #include "AliTPCcalibDB.h"
29 : #include "AliHLTTPCFastTransformObject.h"
30 :
31 : #include <iostream>
32 : #include <iomanip>
33 :
34 : using namespace std;
35 :
36 6 : ClassImp(AliHLTTPCFastTransform); //ROOT macro for the implementation of ROOT specific class methods
37 :
38 :
39 : AliHLTTPCFastTransform::AliHLTTPCFastTransform()
40 : :
41 6 : fMinInitSec(0),
42 6 : fMaxInitSec(fkNSec),
43 6 : fError(),
44 6 : fInitialisationMode(-1),
45 6 : fOrigTransform(0),
46 6 : fLastTimeStamp(-1),
47 6 : fLastTimeBin(600),
48 6 : fTimeBorder1(100),
49 6 : fTimeBorder2(500),
50 6 : fAlignment(NULL)
51 24 : {
52 : // see header file for class documentation
53 : // or
54 : // refer to README to build package
55 : // or
56 : // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
57 876 : for( Int_t i=0; i<fkNSec; i++)
58 87264 : for( Int_t j=0; j<fkNRows; j++ ) fRows[i][j] = NULL;
59 12 : }
60 :
61 : AliHLTTPCFastTransform::~AliHLTTPCFastTransform()
62 24 : {
63 : // see header file for class documentation
64 6 : DeInit();
65 12 : }
66 :
67 :
68 : void AliHLTTPCFastTransform::DeInit()
69 : {
70 : // Deinitialisation
71 :
72 882 : for( Int_t i=0; i<fkNSec; i++){
73 87264 : for( Int_t j=0; j<fkNRows; j++ ){
74 43200 : delete fRows[i][j];
75 43200 : fRows[i][j] = 0;
76 : }
77 : }
78 6 : fOrigTransform = NULL;
79 6 : fLastTimeStamp = -1;
80 6 : delete[] fAlignment;
81 6 : fAlignment = NULL;
82 6 : fError = "";
83 6 : fInitialisationMode = -1;
84 6 : }
85 :
86 :
87 : Int_t AliHLTTPCFastTransform::Init( AliTPCTransform *transform, Long_t TimeStamp )
88 : {
89 : // Initialisation
90 :
91 0 : DeInit();
92 0 : fInitialisationMode = 1;
93 :
94 0 : AliTPCcalibDB* pCalib=AliTPCcalibDB::Instance();
95 0 : if(!pCalib ) return Error( -1, "AliHLTTPCFastTransform::Init: No TPC calibration instance found");
96 :
97 0 : AliTPCParam *tpcParam = pCalib->GetParameters();
98 0 : if( !tpcParam ) return Error( -2, "AliHLTTPCFastTransform::Init: No TPCParam object found");
99 :
100 0 : if( !transform ) transform = pCalib->GetTransform();
101 0 : if( !transform ) return Error( -3, "AliHLTTPCFastTransform::Init: No TPC transformation found");
102 :
103 0 : tpcParam->Update();
104 0 : tpcParam->ReadGeoMatrices();
105 :
106 :
107 : // at the moment initialise all the rows
108 :
109 0 : int nSec = tpcParam->GetNSector();
110 0 : if( nSec>fkNSec ) nSec = fkNSec;
111 :
112 0 : for( Int_t i=0; i<nSec; i++ ){
113 0 : int nRows = tpcParam->GetNRow(i);
114 0 : if( nRows>fkNRows ) nRows = fkNRows;
115 0 : for( int j=0; j<nRows; j++){
116 0 : if( !fRows[i][j] ) fRows[i][j] = new AliHLTTPCFastTransform::AliRowTransform;
117 0 : if( !fRows[i][j] ) return Error( -4, "AliHLTTPCFastTransform::Init: Not enough memory");
118 : }
119 0 : }
120 :
121 0 : const AliTPCRecoParam *rec = transform->GetCurrentRecoParam();
122 :
123 0 : if( !rec ) return Error( -5, "AliHLTTPCFastTransform::Init: No TPC Reco Param set in transformation");
124 :
125 0 : fOrigTransform = transform;
126 :
127 0 : if( rec->GetUseSectorAlignment() && (!pCalib->HasAlignmentOCDB()) ){
128 :
129 0 : fAlignment = new Float_t [fkNSec*21];
130 0 : for( Int_t iSec=0; iSec<fkNSec; iSec++ ){
131 0 : Float_t *t = fAlignment + iSec*21, *r = t+3, *v = t+12;
132 0 : for( int i=0; i<21; i++ ) t[i]=0.f;
133 0 : r[0] = r[4] = r[8] = 1.f;
134 0 : v[0] = v[4] = v[8] = 1.f;
135 : }
136 :
137 0 : for( Int_t iSec=0; iSec<nSec; iSec++ ){
138 0 : Float_t *t = fAlignment + iSec*21, *r = t+3, *v = t+12;
139 0 : TGeoHMatrix *alignment = tpcParam->GetClusterMatrix( iSec );
140 0 : if ( alignment ){
141 0 : const Double_t *tr = alignment->GetTranslation();
142 0 : const Double_t *rot = alignment->GetRotationMatrix();
143 0 : if( tr && rot ){
144 0 : for( int i=0; i<3; i++ ) t[i] = tr[i];
145 0 : for( int i=0; i<9; i++ ) r[i] = rot[i];
146 0 : CalcAdjugateRotation(r,v,1);
147 0 : }
148 0 : }
149 : }
150 0 : }
151 :
152 0 : return SetCurrentTimeStamp( TimeStamp );
153 0 : }
154 :
155 :
156 : Int_t AliHLTTPCFastTransform::WriteToObject( AliHLTTPCFastTransformObject &obj )
157 : {
158 : //
159 : // write fast transformation to ROOT object to store it in database
160 : //
161 0 : obj.Reset();
162 :
163 0 : if( obj.GetNSec() < fkNSec ) return Error( -10, "AliHLTTPCFastTransform::WriteToObject: Wrong N Sectors in object");
164 :
165 0 : obj.SetLastTimeBin( fLastTimeBin );
166 0 : obj.SetTimeSplit1( fTimeBorder1 );
167 0 : obj.SetTimeSplit2( fTimeBorder2 );
168 :
169 0 : TArrayF &alignment =obj.GetAlignmentNonConst();
170 0 : if( !fAlignment ) alignment.Set(0);
171 : else{
172 0 : alignment.Set( fkNSec*21 );
173 0 : for( Int_t i=0; i<fkNSec*21; i++ ) alignment[i] = fAlignment[i];
174 : }
175 :
176 0 : for( Int_t iSec=fMinInitSec; iSec<fMaxInitSec; iSec++ ){
177 0 : for( int iRow=0; iRow<fkNRows; iRow++){
178 0 : if( !fRows[iSec][iRow] ) break;
179 0 : if( ( iSec<obj.GetNSecIn() && iRow>=obj.GetNRowsIn() ) || ( iSec>=obj.GetNSecIn() && iRow>=obj.GetNRowsOut() ) ){
180 0 : return Error( -10, "AliHLTTPCFastTransform::WriteToObject: Wrong N Rows in object");
181 : }
182 0 : for( int iSpline=0; iSpline<3; iSpline++ ){
183 0 : AliHLTTPCSpline2D3D & spline = fRows[iSec][iRow]->fSpline[iSpline];
184 0 : AliHLTTPCSpline2D3DObject& splineObj = obj.GetSplineNonConst( iSec, iRow, iSpline );
185 0 : spline.WriteToObject( splineObj );
186 : }
187 : }
188 0 : obj.SetInitSec(iSec, true);
189 : }
190 :
191 0 : return 0;
192 0 : }
193 :
194 :
195 :
196 : Int_t AliHLTTPCFastTransform::ReadFromObject( const AliHLTTPCFastTransformObject &obj )
197 : {
198 : //
199 : // read fast transformation from ROOT object in database
200 : //
201 :
202 0 : DeInit();
203 0 : fInitialisationMode = 0;
204 :
205 0 : if( obj.GetNSec() > fkNSec ) return Error( -10, "AliHLTTPCFastTransform::ReadFromObject: Wrong N Sectors in object");
206 :
207 0 : fOrigTransform = NULL;
208 0 : fLastTimeStamp = 0;
209 :
210 0 : fLastTimeBin = obj.GetLastTimeBin();
211 0 : fTimeBorder1 = obj.GetTimeSplit1();
212 0 : fTimeBorder2 = obj.GetTimeSplit2();
213 :
214 0 : AliTPCcalibDB* pCalib=AliTPCcalibDB::Instance();
215 0 : if(!pCalib ) return Error( -1, "AliHLTTPCFastTransform::Init: No TPC calibration instance found");
216 :
217 0 : AliTPCParam *tpcParam = pCalib->GetParameters();
218 0 : if( !tpcParam ) return Error( -2, "AliHLTTPCFastTransform::Init: No TPCParam object found");
219 :
220 0 : tpcParam->Update();
221 0 : tpcParam->ReadGeoMatrices();
222 :
223 0 : int nSec = tpcParam->GetNSector();
224 0 : if( nSec>fkNSec ) nSec = fkNSec;
225 :
226 0 : for( Int_t iSec=0; iSec<nSec; iSec++ ){
227 :
228 0 : if (obj.IsSectorInit(iSec) == false) return Error( -10, "AliHLTTPCFastTransform::ReadFromObject: Cannot initialize from Object that does not contain full transform map!");
229 :
230 0 : int nRows = tpcParam->GetNRow(iSec);
231 :
232 0 : if( nRows>fkNRows ) return Error( -10, Form("AliHLTTPCFastTransform::ReadFromObject: N of rows in the fast transformation is too low: %d, needed %d for sector %d", fkNRows, nRows , iSec) );
233 :
234 0 : if( iSec<obj.GetNSecIn() && nRows>obj.GetNRowsIn() ) return Error( -10, Form( "AliHLTTPCFastTransform::ReadFromObject: N of rows in the fast transformation obect is too low: %d, needed %d for sector %d", obj.GetNRowsIn(), nRows , iSec) );
235 :
236 0 : if( iSec>=obj.GetNSecIn() && nRows>obj.GetNRowsOut() ) return Error( -10, Form("AliHLTTPCFastTransform::ReadFromObject: N of rows in the fast transformation obect is too low: %d, needed %d for sector %d", obj.GetNRowsOut(), nRows , iSec) );
237 :
238 0 : for( int iRow=0; iRow<nRows; iRow++){
239 0 : if( !fRows[iSec][iRow] ) fRows[iSec][iRow] = new AliHLTTPCFastTransform::AliRowTransform;
240 0 : if( !fRows[iSec][iRow] ) return Error( -4, "AliHLTTPCFastTransform::ReadFromObject: Not enough memory");
241 0 : for( int iSpline=0; iSpline<3; iSpline++ ){
242 0 : AliHLTTPCSpline2D3D & spline = fRows[iSec][iRow]->fSpline[iSpline];
243 0 : const AliHLTTPCSpline2D3DObject& splineObj = obj.GetSpline( iSec, iRow, iSpline );
244 0 : spline.ReadFromObject( splineObj );
245 : }
246 : }
247 0 : }
248 :
249 0 : const TArrayF &alignment =obj.GetAlignment();
250 0 : delete[] fAlignment;
251 0 : fAlignment = 0;
252 0 : if( alignment.GetSize() == fkNSec*21 ){
253 0 : fAlignment = new Float_t [fkNSec*21];
254 0 : for( Int_t i=0; i<fkNSec*21; i++ ) fAlignment[i] = alignment[i];
255 0 : }
256 : return 0;
257 0 : }
258 :
259 :
260 :
261 : bool AliHLTTPCFastTransform::CalcAdjugateRotation(const Float_t *mA, Float_t *mB, bool bCheck)
262 : {
263 : // check rotation matrix and adjugate for consistency
264 : //
265 : // ( for a rotation matrix inverse == transpose )
266 : //
267 :
268 0 : mB[0] = mA[0];
269 0 : mB[1] = mA[3];
270 0 : mB[2] = mA[6];
271 :
272 0 : mB[3] = mA[1];
273 0 : mB[4] = mA[4];
274 0 : mB[5] = mA[7];
275 :
276 0 : mB[6] = mA[2];
277 0 : mB[7] = mA[5];
278 0 : mB[8] = mA[8];
279 :
280 0 : if (bCheck) {
281 0 : for (int r=0; r<3; r++) {
282 0 : for (int c=0; c<3; c++) {
283 : float a=0.;
284 : float expected=0.;
285 0 : if (r==c) expected=1.;
286 0 : for (int i=0; i<3; i++) {
287 0 : a+=mA[3*r+i]*mB[c+(3*i)];
288 : }
289 0 : if (TMath::Abs(a-expected)>0.00001) {
290 0 : std::cout << "inconsistent adjugate at " << r << c << ": " << a << " " << expected << std::endl;
291 0 : for( int i=0; i<9; i++ ) mB[i] = 0;
292 0 : mB[0] = mB[4] = mB[8] = 1;
293 0 : return false;
294 : }
295 0 : }
296 : }
297 : }
298 0 : return true;
299 0 : }
300 :
301 :
302 :
303 : Int_t AliHLTTPCFastTransform::SetCurrentTimeStamp( Long_t TimeStamp )
304 : {
305 : // Set the current time stamp
306 :
307 0 : if( fInitialisationMode!=1 ){
308 0 : fLastTimeStamp = TimeStamp;
309 0 : return 0;
310 : }
311 :
312 0 : Long_t lastTS = fLastTimeStamp;
313 0 : fLastTimeStamp = -1; // deinitialise
314 :
315 0 : if( !fOrigTransform ) return Error( -1, "AliHLTTPCFastTransform::SetCurrentTimeStamp: TPC transformation has not been set properly");
316 :
317 0 : AliTPCcalibDB* pCalib=AliTPCcalibDB::Instance();
318 0 : if(!pCalib ) return Error( -2, "AliHLTTPCFastTransform::SetCurrentTimeStamp: No TPC calibration found");
319 :
320 0 : AliTPCParam *tpcParam = pCalib->GetParameters();
321 0 : if( !tpcParam ) return Error( -3, "AliHLTTPCFastTransform::SetCurrentTimeStamp: No TPCParam object found");
322 :
323 :
324 0 : if( TimeStamp<0 ) return 0;
325 :
326 0 : fLastTimeStamp = lastTS;
327 :
328 0 : if( fLastTimeStamp>=0 && TMath::Abs(fLastTimeStamp - TimeStamp ) <60 ) return 0;
329 :
330 :
331 0 : fOrigTransform->SetCurrentTimeStamp( static_cast<UInt_t>(TimeStamp) );
332 0 : fLastTimeStamp = TimeStamp;
333 :
334 : // find last calibrated time bin
335 :
336 0 : Int_t nTimeBins = tpcParam->GetMaxTBin();
337 0 : Int_t is[]={0};
338 : bool sign = 0;
339 0 : for( fLastTimeBin=0; fLastTimeBin<nTimeBins; fLastTimeBin++){
340 : // static cast is okay since fLastTimeBin has limited value range
341 0 : Double_t xx[]={0,0,static_cast<Double_t>(fLastTimeBin)};
342 0 : fOrigTransform->Transform(xx,is,0,1);
343 0 : bool s = (xx[2]>=0);
344 0 : if( fLastTimeBin==0 ) sign = s;
345 0 : else if( sign!=s ){
346 0 : fLastTimeBin--;
347 0 : break;
348 : }
349 0 : }
350 0 : fTimeBorder1 = 60;
351 0 : fTimeBorder2 = fLastTimeBin - 100;
352 :
353 0 : int nSec = tpcParam->GetNSector();
354 0 : if( nSec>fkNSec ) nSec = fkNSec;
355 :
356 0 : for( Int_t i=fMinInitSec; i<fMaxInitSec; i++ ){
357 0 : int nRows = tpcParam->GetNRow(i);
358 0 : if( nRows>fkNRows ) nRows = fkNRows;
359 0 : for( int j=0; j<nRows; j++){
360 0 : if( fRows[i][j] ){
361 0 : int err = InitRow(i,j);
362 0 : if( err!=0 ) return err;
363 0 : }
364 : }
365 0 : }
366 0 : return 0;
367 0 : }
368 :
369 :
370 : Int_t AliHLTTPCFastTransform::InitRow( Int_t iSector, Int_t iRow )
371 : {
372 : // see header file for class documentation
373 :
374 0 : AliTPCcalibDB* pCalib=AliTPCcalibDB::Instance();
375 :
376 0 : if( iSector<0 || iSector>=fkNSec || iRow<0 || iRow>=fkNRows || !fOrigTransform || (fLastTimeStamp<0) ||
377 0 : !fRows[iSector][iRow] || !pCalib || !pCalib->GetParameters() ){
378 0 : return Error( -1, "AliHLTTPCFastTransform::InitRow: Internal error");
379 : }
380 :
381 0 : AliTPCParam *tpcParam = pCalib->GetParameters();
382 :
383 0 : Int_t nPads = tpcParam->GetNPads(iSector,iRow);
384 :
385 0 : if( nPads<2 ){
386 0 : return Error( -2, Form("AliHLTTPCFastTransform::InitRow: Wrong NPads=%d for sector %d row %d",nPads,iSector,iRow ));
387 : }
388 :
389 0 : fRows[iSector][iRow]->fSpline[0].Init( 0.5, nPads-1+0.5, 15, 0, fTimeBorder1, 5);
390 0 : fRows[iSector][iRow]->fSpline[1].Init( 0.5, nPads-1+0.5, 15, fTimeBorder1, fTimeBorder2, 10);
391 0 : fRows[iSector][iRow]->fSpline[2].Init( 0.5, nPads-1+0.5, 15, fTimeBorder2, fLastTimeBin, 5);
392 :
393 0 : for( Int_t i=0; i<3; i++){
394 0 : Int_t is[]={iSector};
395 0 : for( Int_t j=0; j<fRows[iSector][iRow]->fSpline[i].GetNPoints(); j++){
396 0 : Float_t pad, time;
397 0 : fRows[iSector][iRow]->fSpline[i].GetAB(j,pad,time);
398 0 : Double_t xx[]={static_cast<Double_t>(iRow),pad,time};
399 0 : fOrigTransform->Transform(xx,is,0,1);
400 0 : fRows[iSector][iRow]->fSpline[i].Fill(j,xx);
401 0 : }
402 0 : fRows[iSector][iRow]->fSpline[i].Consolidate();
403 0 : }
404 0 : return 0;
405 0 : }
406 :
407 :
408 :
409 : Int_t AliHLTTPCFastTransform::GetRowSize( Int_t iSec, Int_t iRow ) const
410 : {
411 : // see header file for class documentation
412 : Int_t s = sizeof(AliHLTTPCFastTransform::AliRowTransform);
413 0 : if( fRows[iSec][iRow] ) for( Int_t i=0; i<3; i++) s+=fRows[iSec][iRow]->fSpline[i].GetMapSize();
414 0 : return s;
415 : }
416 :
417 : Int_t AliHLTTPCFastTransform::GetSize() const
418 : {
419 : // see header file for class documentation
420 : Int_t s = sizeof(AliHLTTPCFastTransform);
421 0 : for( Int_t i=0; i<fkNSec; i++ )
422 0 : for( Int_t j=0; j<fkNRows; j++ ) if( fRows[i][j] ){
423 0 : s+= sizeof(AliHLTTPCFastTransform::AliRowTransform);
424 0 : for( Int_t k=0; k<3; k++) fRows[i][j]->fSpline[k].GetMapSize();
425 0 : }
426 0 : return s;
427 : }
428 :
429 :
430 : void AliHLTTPCFastTransform::Print(const char* /*option*/) const
431 : {
432 : // print info
433 0 : ios::fmtflags coutflags=std::cout.flags(); // backup cout status flags
434 :
435 0 : if( !fAlignment ){
436 0 : std::cout << "AliHLTTPCFastTransform: no alignment transformation";
437 0 : } else {
438 0 : for( int iSec=0; iSec<fkNSec; iSec++ ){
439 0 : std::cout << "AliHLTTPCClusterTransformation for sector " << iSec << std::endl;
440 :
441 0 : const Float_t *mT = fAlignment + iSec*21;
442 0 : const Float_t *mR = mT + 3;
443 0 : const Float_t *mV = mT + 3+9;
444 :
445 0 : std::cout.setf(ios_base::showpos|ios_base::showpos|ios::right);
446 0 : std::cout << " translation: " << std::endl;
447 : int r=0;
448 0 : for (r=0; r<3; r++) {
449 0 : std::cout << setw(7) << fixed << setprecision(2);
450 0 : cout << " " << mT[r] << std::endl;
451 : }
452 0 : std::cout << " rotation and adjugated rotation: " << std::endl;
453 0 : for (r=0; r<3; r++) {
454 : int c=0;
455 0 : std::cout << setw(7) << fixed << setprecision(2);
456 0 : for (c=0; c<3; c++) std::cout << " " << mR[3*r+c];
457 0 : std::cout << " ";
458 0 : for (c=0; c<3; c++) std::cout << " " << mV[3*r+c];
459 0 : std::cout << endl;
460 : }
461 : }
462 : }
463 0 : std::cout.flags(coutflags); // restore the original flags
464 0 : }
|
