15 minutes interval for calculation of drift correction
[u/mrichter/AliRoot.git] / TPC / AliTPCTempMap.cxx
CommitLineData
1209231c 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16
17///////////////////////////////////////////////////////////////////////////////
18// //
19// TPC calibration class for temperature maps and tendencies //
20// (based on TPC Temperature Sensors and FiniteElement Simulation) //
21// //
22// Authors: Stefan Rossegger, Haavard Helstrup //
23// //
24///////////////////////////////////////////////////////////////////////////////
25
26#include "AliTPCSensorTempArray.h"
27#include "TLinearFitter.h"
28#include "TString.h"
29#include "TGraph2D.h"
f59cd9d0 30#include "TTimeStamp.h"
1209231c 31
32#include "AliTPCTempMap.h"
33
34
35ClassImp(AliTPCTempMap)
36
37 const char kStringFEsimulation[] = "FEsimulation.txt";
38
39//_____________________________________________________________________________
40AliTPCTempMap::AliTPCTempMap(AliTPCSensorTempArray *sensorDCS):
41 TNamed(),
f59cd9d0 42 fTempArray(0),
1209231c 43 fStringFEsimulation(kStringFEsimulation)
44{
45 //
46 // AliTPCTempMap default constructor
47 //
48
f59cd9d0 49 fTempArray = sensorDCS;
1209231c 50
51}
52
53//_____________________________________________________________________________
54AliTPCTempMap::AliTPCTempMap(const AliTPCTempMap &c):
55 TNamed(c),
f59cd9d0 56 fTempArray(c.fTempArray),
1209231c 57 fStringFEsimulation(c.fStringFEsimulation)
58{
59 //
60 // AliTPCTempMap copy constructor
61 //
62
63}
64
65//_____________________________________________________________________________
66AliTPCTempMap::~AliTPCTempMap()
67{
68 //
69 // AliTPCTempMap destructor
70 //
71
72}
73
74//_____________________________________________________________________________
75AliTPCTempMap &AliTPCTempMap::operator=(const AliTPCTempMap &c)
76{
77 //
78 // Assignment operator
79 //
80
81 if (this != &c) ((AliTPCTempMap &) c).Copy(*this);
82 return *this;
83
84}
85
86//_____________________________________________________________________________
87void AliTPCTempMap::Copy(TObject &c) const
88{
89 //
90 // Copy function
91 //
92
93 TObject::Copy(c);
94
95}
96
97//_____________________________________________________________________________
98
99Double_t AliTPCTempMap::GetTempGradientY(UInt_t timeSec, Int_t side){
100 //
101 // Extract Linear Vertical Temperature Gradient [K/cm] within the TPC on
102 // Shaft Side(A): 0
103 // Muon Side(C): 1
104 // Values based on TemperatureSensors within the TPC ( type: 3 (TPC) )
105 //
106 // FIXME: Also return residual-distribution, covariance Matrix
107 // or simply chi2 for validity check?
108 // -> better use GetLinearFitter - function in this case!
109
b479e253 110 TLinearFitter *fitter = new TLinearFitter(3,"x0++x1++x2");
1209231c 111 TVectorD param(3);
112 Int_t i = 0;
113
f59cd9d0 114 Int_t nsensors = fTempArray->NumSensors();
1209231c 115 for (Int_t isensor=0; isensor<nsensors; isensor++) { // loop over all sensors
f59cd9d0 116 AliTPCSensorTemp *entry = (AliTPCSensorTemp*)fTempArray->GetSensorNum(isensor);
1209231c 117
118 if (entry->GetType()==3 && entry->GetSide()==side) { // take SensorType:TPC
119 Double_t x[3];
120 x[0]=1;
121 x[1]=entry->GetX();
122 x[2]=entry->GetY();
f59cd9d0 123 Double_t y = fTempArray->GetValue(timeSec,isensor); // get temperature value
f4fe7830 124 if (IsOK(y)) fitter->AddPoint(x,y,1); // add values to LinearFitter
1209231c 125 i++;
126 }
127
128 }
b479e253 129 fitter->Eval();
130 fitter->GetParameters(param);
131
132 fitter->~TLinearFitter();
1209231c 133
134 return param[2]; // return vertical (Y) tempGradient in [K/cm]
135
136}
137
f1ea1647 138//_____________________________________________________________________________
f59cd9d0 139TLinearFitter *AliTPCTempMap::GetLinearFitter(Int_t type, Int_t side, TTimeStamp &stamp)
140{
141 //
142 // absolute time stamp used
143 // see AliTPCTempMap::GetLinearFitter(Int_t type, Int_t side, UInt_t timeSec) for details
144 //
148518a2 145 Int_t timeSec = stamp.GetSec()-fTempArray->GetStartTime().GetSec();
f59cd9d0 146 return GetLinearFitter(type,side,timeSec);
147}
148
1209231c 149//_____________________________________________________________________________
1209231c 150TLinearFitter *AliTPCTempMap::GetLinearFitter(Int_t type, Int_t side, UInt_t timeSec)
151{
152 //
153 // Creates a TlinearFitter object for the desired region of the TPC
154 // (via choosen type and side of TPC temperature sensors) at a given
155 // timeSec (in secounds) after start time
156 // type: 0 ... ReadOutChambers (ROC)
157 // 1 ... OuterContainmentVessel (OFC)
158 // 2 ... InnerContainmentVessel (IFC) + ThermalScreener (TS)
159 // 3 ... Within the TPC (DriftVolume) (TPC)
d0bd4fcc 160 // 4 ... Only InnerContainmentVessel (IFC)
1209231c 161 // side: Can be choosen for type 0 and 3 (otherwise it will be ignored in
162 // in order to get all temperature sensors of interest)
163 // 0 ... Shaft Side (A)
164 // 1 ... Muon Side (C)
165 //
166
167 TLinearFitter *fitter = new TLinearFitter(3);
168 Double_t *x = new Double_t[3];
169 Double_t y = 0;
f4fe7830 170 const Float_t kMaxDelta=0.5;
f59cd9d0 171
d0bd4fcc 172 if (type == 1 || type == 2 || type == 4) {
1209231c 173 fitter->SetFormula("x0++x1++TMath::Sin(x2)"); // returns Z,Y gradient
174 } else {
175 fitter->SetFormula("x0++x1++x2"); // returns X,Y gradient
176 }
177
178 Int_t i = 0;
f59cd9d0 179 Int_t nsensors = fTempArray->NumSensors();
180
181 Float_t temps[1000];
182 for (Int_t isensor=0; isensor<nsensors; isensor++) { // loop over all sensors
183 AliTPCSensorTemp *entry = (AliTPCSensorTemp*)fTempArray->GetSensorNum(isensor);
184 if (entry->GetType()==type && entry->GetSide()==side){
185 Float_t temperature= fTempArray->GetValue(timeSec,isensor); // get temperature value
f4fe7830 186 if (IsOK(temperature)) {temps[i]=temperature; i++;}
f59cd9d0 187 }
188 }
189 Float_t medianTemp = TMath::Median(i, temps);
190 if (i<3) return 0;
191 Float_t rmsTemp = TMath::RMS(i, temps);
192
193 i=0;
194
1209231c 195 for (Int_t isensor=0; isensor<nsensors; isensor++) { // loop over all sensors
f59cd9d0 196 AliTPCSensorTemp *entry = (AliTPCSensorTemp*)fTempArray->GetSensorNum(isensor);
1209231c 197
198 if (type==0 || type==3) { // 'side' information used
199 if (entry->GetType()==type && entry->GetSide()==side) {
200 x[0]=1;
201 x[1]=entry->GetX();
202 x[2]=entry->GetY();
f59cd9d0 203 y = fTempArray->GetValue(timeSec,isensor); // get temperature value
204 if (TMath::Abs(y-medianTemp)>kMaxDelta+4.*rmsTemp) continue;
f4fe7830 205 if (IsOK(y)) fitter->AddPoint(x,y,1); // add values to LinearFitter
1209231c 206 i++;
207 }
208 } else if (type==2) { // in case of IFC also usage of TS values
209 if ((entry->GetType()==2) || (entry->GetType()==5)) {
210 x[0]=1;
211 x[1]=entry->GetZ();
212 x[2]=entry->GetPhi();
f59cd9d0 213 y = fTempArray->GetValue(timeSec,isensor);
214 if (TMath::Abs(y-medianTemp)>kMaxDelta+4.*rmsTemp) continue;
f4fe7830 215 if (IsOK(y)) fitter->AddPoint(x,y,1);
1209231c 216 i++;
217 }
218 } else if (type==1){
219 if (entry->GetType()==type) {
220 x[0]=1;
221 x[1]=entry->GetZ();
222 x[2]=entry->GetPhi();
f59cd9d0 223 y = fTempArray->GetValue(timeSec,isensor);
224 if (TMath::Abs(y-medianTemp)>kMaxDelta+4.*rmsTemp) continue;
f4fe7830 225 if (IsOK(y)) fitter->AddPoint(x,y,1);
1209231c 226 i++;
227 }
d0bd4fcc 228 } else if (type==4) { // ONLY IFC
229 if (entry->GetType()==2) {
230 x[0]=1;
231 x[1]=entry->GetZ();
232 x[2]=entry->GetPhi();
f59cd9d0 233 y = fTempArray->GetValue(timeSec,isensor);
234 if (TMath::Abs(y-medianTemp)>kMaxDelta+4.*rmsTemp) continue;
f4fe7830 235 if (IsOK(y)) fitter->AddPoint(x,y,1);
d0bd4fcc 236 i++;
237 }
1209231c 238 }
239 }
f59cd9d0 240 fitter->Eval();
241 //fitter->EvalRobust(0.9); // Evaluates fitter
1209231c 242
243 delete [] x;
244
245 return fitter;
246
247 // returns TLinearFitter object where Chi2, Fitparameters and residuals can
248 // be extracted via usual memberfunctions
b479e253 249 // example: fitter->GetParameters(param)
1209231c 250 // In case of type IFC or OFC, the parameters are the gradients in
251 // Z and Y direction (see fitformula)
252 // Caution: Parameters are [K/cm] except Y at IFC,OFC ([K/radius])
253}
254
255//_____________________________________________________________________________
256
257TGraph2D *AliTPCTempMap::GetTempMapsViaSensors(Int_t type, Int_t side, UInt_t timeSec)
258{
259 //
260 // Creates a TGraph2D object for the desired region of the TPC
261 // (via choosen type and side of TPC temperature sensors) at a given
262 // timeSec (in secounds) after start time
263 // type: 0 ... ReadOutChambers (ROC)
264 // 1 ... OuterContainmentVessel (OFC)
265 // 2 ... InnerContainmentVessel (IFC) + ThermalScreener (TS)
266 // 3 ... Within the TPC (DriftVolume) (TPC)
267 // side: Can be choosen for type 0 and 3 (otherwise it will be ignored in
268 // in order to get all temperature sensors of interest)
55f06b51 269 // 0 ... A - side
270 // 1 ... C - side
1209231c 271 //
272
273 TGraph2D *graph2D = new TGraph2D();
274
275 Int_t i = 0;
276
f59cd9d0 277 Int_t nsensors = fTempArray->NumSensors();
1209231c 278 for (Int_t isensor=0; isensor<nsensors; isensor++) { // loop over all sensors
f59cd9d0 279 AliTPCSensorTemp *entry = (AliTPCSensorTemp*)fTempArray->GetSensorNum(isensor);
1209231c 280
281 Double_t x, y, z, r, phi, tempValue;
282 x = entry->GetX();
283 y = entry->GetY();
284 z = entry->GetZ();
285 r = entry->GetR();
286 phi = entry->GetPhi();
f59cd9d0 287 tempValue = fTempArray->GetValue(timeSec,isensor);
55f06b51 288 // printf("%d type %d: x=%lf y=%lf temp=%lf\n",isensor,entry->GetType(),x,y, tempValue);
1209231c 289 if (type==0 || type==3) { // 'side' information used
290 if (entry->GetType()==type && entry->GetSide()==side) {
291 graph2D->SetPoint(i,x,y,tempValue);
292 i++;
293 }
294 } else if (type==2) { // in case of IFC also usage of TS values
295 if (entry->GetType()==2 || entry->GetType()==5) {
296 graph2D->SetPoint(i,z,phi,tempValue);
297 i++;
298 }
299 } else if (type==1){
300 if (entry->GetType()==type) {
301 graph2D->SetPoint(i,z,phi,tempValue);
302 i++;
303 }
304 }
305 }
306
307 if (type==0 || type==3) {
308 graph2D->GetXaxis()->SetTitle("X[cm]");
309 graph2D->GetYaxis()->SetTitle("Y[cm]");
310 if (type==0 && side==0) {
55f06b51 311 graph2D->SetTitle("ROC A side");
1209231c 312 } else if (type==0 && side==1) {
55f06b51 313 graph2D->SetTitle("ROC C side");
1209231c 314 } else if (type==3 && side==0) {
55f06b51 315 graph2D->SetTitle("TPC A side (Inside the TPC)");
1209231c 316 } else if (type==3 && side==1) {
55f06b51 317 graph2D->SetTitle("TPC C side (Inside the TPC)");
1209231c 318 }
319 } else if (type==1 || type==2) {
320 graph2D->GetXaxis()->SetTitle("Z[cm]");
321 graph2D->GetYaxis()->SetTitle("Phi[RAD]");
322 if (type==1) {
323 graph2D->SetTitle("Outer Containment Vessel");
324 } else if (type==2) {
55f06b51 325 graph2D->SetTitle("Inner Containment Vessel");
1209231c 326 }
327 }
328
329 if (!graph2D->GetN()) {
330 printf("Returned TGraph2D is empty: check type and side values\n");
331 }
332
333 graph2D->GetXaxis()->SetLabelOffset(0.0);
334 graph2D->GetYaxis()->SetLabelOffset(0.005);
335 graph2D->GetZaxis()->SetLabelOffset(-0.04);
336
337
338 return graph2D; // returns TGgraph2D object
339
340}
341
342
343//_____________________________________________________________________________
344
345TGraph *AliTPCTempMap::MakeGraphGradient(Int_t axis, Int_t side, Int_t nPoints)
346{
347 //
348 // Make graph from start time to end time of TempGradient in axis direction
349 // axis: 0 ... horizontal Temperature Gradient (X)
350 // 1 ... vertical Temperature Gradient (Y)
351 // 2 ... longitudenal Temperature Gradient (Z) (side is ignored)
352 // z gradient value based on OFC temperature sensors
353 // Caution!: better z gradient values through difference between
354 // param[0] A- and param[0] C-side !
355 // side for X and Y gradient:
356 // 0 ... Shaft Side (A)
357 // 1 ... Muon Side (C)
358 //
359
360 TVectorD param(3);
361 TLinearFitter *fitter = new TLinearFitter(3);
362
f59cd9d0 363 UInt_t fStartTime = fTempArray->AliTPCSensorTempArray::GetStartTime();
364 UInt_t fEndTime = fTempArray->AliTPCSensorTempArray::GetEndTime();
1209231c 365
366 UInt_t stepTime = (fEndTime-fStartTime)/nPoints;
367
368 Double_t *x = new Double_t[nPoints];
369 Double_t *y = new Double_t[nPoints];
370 for (Int_t ip=0; ip<nPoints; ip++) {
371 x[ip] = fStartTime+ip*stepTime;
372 if (axis==2) {// Gradient in Z direction (based on OFC tempSensors)
373 fitter = GetLinearFitter(1, side, ip*stepTime);
374 } else {// Gradient in X or Y direction (based on TPC tempSensors)
375 fitter = GetLinearFitter(3, side, ip*stepTime);
376 }
377 fitter->GetParameters(param);
378 // multiplied by 500 since TempGradient is in [K/cm]
379 // (TPC diameter and length ~500cm)
380 if (axis==1) { // Y axis
381 y[ip] = param[2]*500;
382 } else { // X axis
383 y[ip] = param[1]*500;
384 }
385 }
386
387 TGraph *graph = new TGraph(nPoints,x,y);
388
389 fitter->~TLinearFitter();
390 delete [] x;
391 delete [] y;
392
393 graph->GetXaxis()->SetTimeDisplay(1);
394 graph->GetXaxis()->SetLabelOffset(0.02);
395 graph->GetXaxis()->SetTimeFormat("#splitline{%d/%m}{%H:%M}");
396
397 return graph;
398}
399
f1ea1647 400
401//_____________________________________________________________________________
402Double_t AliTPCTempMap::GetTemperature(Double_t x, Double_t y, Double_t z, TTimeStamp &stamp)
403{
404 //
405 // absolute time stamp used
406 // see also Double_t AliTPCTempMap::GetTemperature(Double_t x, Double_t y, Double_t z, UInt_t timeSec) for details
407 //
408
409 Int_t timeSec = stamp.GetSec()-fTempArray->GetStartTime().GetSec();
410 return GetTemperature(x, y, z, timeSec);
411}
412
1209231c 413//_____________________________________________________________________________
414
415Double_t AliTPCTempMap::GetTemperature(Double_t x, Double_t y, Double_t z, UInt_t timeSec)
416{
417 //
f1ea1647 418 // Returns estimated Temperature at given position (x,y,z[cm]) at given time
1209231c 419 // (timeSec) after starttime
420 // Method: so far just a linear interpolation between Linar fits of
421 // the TPC temperature sensors
422 // FIXME: 'Educated Fit' through FiniteElement Simulation results!
423 // FIXXME: Return 0? if x,y,z out of range
424 //
425
426 TVectorD paramA(3), paramC(3);
427 TLinearFitter *fitterA = new TLinearFitter(3);
428 TLinearFitter *fitterC = new TLinearFitter(3);
429
430 fitterA = GetLinearFitter(3, 0, timeSec);
431 fitterA->GetParameters(paramA);
432 fitterC = GetLinearFitter(3, 1, timeSec);
433 fitterC->GetParameters(paramC);
434
435 Double_t fvalA = paramA[0]+paramA[1]*x+paramA[2]*y;
436 Double_t fvalC = paramC[0]+paramC[1]*x+paramC[2]*y;
437
438 Double_t k = (fvalA-fvalC)/(2*247);
439 Double_t tempValue = fvalC+(fvalA-fvalC)/2+k*z;
440
441 fitterA->~TLinearFitter();
442 fitterC->~TLinearFitter();
443
444 return tempValue;
f1ea1647 445
1209231c 446}
447
f4fe7830 448
449Bool_t AliTPCTempMap::IsOK(Float_t value){
450 //
451 //
452 //
453 const Float_t kMinT=15;
454 const Float_t kMaxT=25;
455 return (value>kMinT && value<kMaxT);
456}