1 /**************************************************************************
2 * Copyright(c) 1998-2011, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
19 // Mean Vertex preprocessor:
20 // 2) takes data after pass0 ,
21 // processes it, and stores either to OCDB .
25 #include "AliMeanVertexPreprocessorOffline.h"
27 #include "AliCDBStorage.h"
28 #include "AliCDBMetaData.h"
29 #include "AliCDBManager.h"
30 #include "AliCDBEntry.h"
33 #include <TTimeStamp.h>
35 #include <TObjString.h>
40 #include "AliESDVertex.h"
47 ClassImp(AliMeanVertexPreprocessorOffline)
49 //_______________________________________________________
51 const Char_t *AliMeanVertexPreprocessorOffline::fgkStatusCodeName[AliMeanVertexPreprocessorOffline::kNStatusCodes] = {
53 "open file error or missing histos",
55 "problems storing OCDB",
56 "write MeanVertex computed online",
57 "write SPD vtx offline",
58 "lumi region or cov matrix computation problems, default values set"
62 //____________________________________________________
63 AliMeanVertexPreprocessorOffline::AliMeanVertexPreprocessorOffline():
64 TNamed("AliMeanVertexPreprocessorOffline","AliMeanVertexPreprocessorOffline"),
69 //____________________________________________________
71 AliMeanVertexPreprocessorOffline::~AliMeanVertexPreprocessorOffline()
76 //____________________________________________________
77 void AliMeanVertexPreprocessorOffline::ProcessOutput(const char *filename, AliCDBStorage *db, Int_t runNb){
79 TFile *file = TFile::Open(filename);
80 if (!file || !file->IsOpen()){
81 AliError(Form("cannot open outputfile %s", filename));
87 AliError(Form("no OCDB storage found, return"));
93 TList *list = (TList*)file->Get("MeanVertex");
95 TH1F *histTRKvtxX = 0x0;
96 TH1F *histTRKvtxY = 0x0;
97 TH1F *histTRKvtxZ = 0x0;
99 TH1F *histSPDvtxX = 0x0;
100 TH1F *histSPDvtxY = 0x0;
101 TH1F *histSPDvtxZ = 0x0;
104 Bool_t useTRKvtx = kTRUE;
105 Bool_t useITSSAvtx = kFALSE;
106 Bool_t useSPDvtx = kFALSE;
107 Bool_t spdAvailable = kTRUE;
108 Bool_t writeMeanVertexSPD = kFALSE;
109 Bool_t vertexerSPD3Doff=kFALSE;
114 histTRKvtxX = (TH1F*)file->Get("hTRKVertexX");
115 histTRKvtxY = (TH1F*)file->Get("hTRKVertexY");
116 histTRKvtxZ = (TH1F*)file->Get("hTRKVertexZ");
118 histSPDvtxX = (TH1F*)file->Get("hSPDVertexX");
119 histSPDvtxY = (TH1F*)file->Get("hSPDVertexY");
120 histSPDvtxZ = (TH1F*)file->Get("hSPDVertexZ");
122 if (!histTRKvtxX || !histTRKvtxY || !histTRKvtxZ){
127 histTRKvtxX = (TH1F*)file->FindObject("hITSSAVertexX");
128 histTRKvtxY = (TH1F*)file->FindObject("hITSSAVertexY");
129 histTRKvtxZ = (TH1F*)file->FindObject("hITSSAVertexZ");
131 if (!histTRKvtxX || !histTRKvtxY || !histTRKvtxZ){
136 if (!histSPDvtxX || !histSPDvtxY || !histSPDvtxZ){
137 AliError(Form("cannot find any histograms available from file"));
147 histTRKvtxX = (TH1F*)list->FindObject("hTRKVertexX");
148 histTRKvtxY = (TH1F*)list->FindObject("hTRKVertexY");
149 histTRKvtxZ = (TH1F*)list->FindObject("hTRKVertexZ");
151 histSPDvtxX = (TH1F*)list->FindObject("hSPDVertexX");
152 histSPDvtxY = (TH1F*)list->FindObject("hSPDVertexY");
153 histSPDvtxZ = (TH1F*)list->FindObject("hSPDVertexZ");
155 if (!histTRKvtxX || !histTRKvtxY || !histTRKvtxZ){
160 histTRKvtxX = (TH1F*)list->FindObject("hITSSAVertexX");
161 histTRKvtxY = (TH1F*)list->FindObject("hITSSAVertexY");
162 histTRKvtxZ = (TH1F*)list->FindObject("hITSSAVertexZ");
164 if (!histTRKvtxX || !histTRKvtxY || !histTRKvtxZ){
169 if (!histSPDvtxX || !histSPDvtxY || !histSPDvtxZ){
170 AliError(Form("cannot find any histograms available from list"));
182 Float_t nEntriesX = histTRKvtxX->GetEffectiveEntries();
183 Float_t nEntriesY = histTRKvtxY->GetEffectiveEntries();
184 Float_t nEntriesZ = histTRKvtxZ->GetEffectiveEntries();
186 if (nEntriesX < 50. || nEntriesY<50. || nEntriesZ<50.) {
187 AliError(Form("TRK vertex histograms have too few entries for fitting"));
194 Float_t nEntriesX = histTRKvtxX->GetEffectiveEntries();
195 Float_t nEntriesY = histTRKvtxY->GetEffectiveEntries();
196 Float_t nEntriesZ = histTRKvtxZ->GetEffectiveEntries();
198 if (nEntriesX < 50. || nEntriesY<50. || nEntriesZ<50.) {
199 AliError(Form("ITSSA vertex histograms have too few entries for fitting"));
205 Float_t nEntriesX = histSPDvtxX->GetEffectiveEntries();
206 Float_t nEntriesY = histSPDvtxY->GetEffectiveEntries();
207 Float_t nEntriesZ = histSPDvtxZ->GetEffectiveEntries();
209 if (nEntriesX < 50. || nEntriesY<50. || nEntriesZ<50.) {
210 spdAvailable = kFALSE;
211 if ((useTRKvtx==kFALSE) && (useITSSAvtx==kFALSE)){
212 AliError(Form("Also SPD vertex histograms have too few entries for fitting, return"));
213 fStatus=kLowStatistics;
218 if((nEntriesX == 0.)&&(nEntriesY==0.) && (nEntriesZ>0.)){
219 vertexerSPD3Doff = kTRUE;
220 AliWarning("Vertexer SPD 3D off");
223 Double_t xMeanVtx=0., yMeanVtx=0., zMeanVtx=0.;
224 Double_t xSigmaVtx=0., ySigmaVtx=0., zSigmaVtx=0.;
227 TF1 *fitVtxX, *fitVtxY, *fitVtxZ;
229 if (useTRKvtx || useITSSAvtx){
230 histTRKvtxX ->Fit("gaus", "M");
231 fitVtxX = histTRKvtxX -> GetFunction("gaus");
232 xMeanVtx = fitVtxX -> GetParameter(1);
233 if (TMath::Abs(xMeanVtx) > 2.) {
235 writeMeanVertexSPD=kTRUE;
236 fStatus=kWriteMeanVertexSPD;
239 histTRKvtxY ->Fit("gaus", "M");
240 fitVtxY = histTRKvtxY -> GetFunction("gaus");
241 yMeanVtx = fitVtxY -> GetParameter(1);
242 if (TMath::Abs(yMeanVtx) > 2.) {
244 writeMeanVertexSPD=kTRUE;
245 fStatus=kWriteMeanVertexSPD;
248 histTRKvtxZ ->Fit("gaus", "M", "", -12, 12);
249 fitVtxZ = histTRKvtxZ -> GetFunction("gaus");
250 zMeanVtx = fitVtxZ -> GetParameter(1);
251 zSigmaVtx = fitVtxZ -> GetParameter(2);
252 if ((TMath::Abs(zMeanVtx) > 20.) || (zSigmaVtx>12.)) {
253 zMeanVtx = histTRKvtxZ->GetMean();
254 zSigmaVtx = histTRKvtxZ->GetRMS();
255 writeMeanVertexSPD=kTRUE;
256 fStatus=kWriteMeanVertexSPD;
262 //check fits: compare histo mean with fit mean value
263 Double_t xHistoMean, yHistoMean, zHistoMean;
264 Double_t xHistoRMS, yHistoRMS, zHistoRMS;
266 if (useTRKvtx || useITSSAvtx){
267 xHistoMean = histTRKvtxX -> GetMean();
268 xHistoRMS = histTRKvtxX ->GetRMS();
270 if ((TMath::Abs(xHistoMean-xMeanVtx) > 0.5)){
271 AliWarning(Form("Possible problems with the fit mean very different from histo mean... using SPD vertex"));
273 useITSSAvtx = kFALSE;
275 fStatus=kUseOfflineSPDvtx;
278 yHistoMean = histTRKvtxY ->GetMean();
279 yHistoRMS = histTRKvtxY ->GetRMS();
281 if ((TMath::Abs(yHistoMean-yMeanVtx) > 0.5)){
282 AliWarning(Form("Possible problems with the fit mean very different from histo mean... using SPD vertex"));
284 useITSSAvtx = kFALSE;
286 fStatus=kUseOfflineSPDvtx;
289 zHistoMean = histTRKvtxZ -> GetMean();
290 zHistoRMS = histTRKvtxZ ->GetRMS();
292 if ((TMath::Abs(zHistoMean-zMeanVtx) > 1.)){
293 AliWarning(Form("Possible problems with the fit mean very different from histo mean... using SPD vertex"));
295 useITSSAvtx = kFALSE;
297 fStatus=kUseOfflineSPDvtx;
302 if ((useSPDvtx) && (spdAvailable) && (!vertexerSPD3Doff)){
304 histSPDvtxX ->Fit("gaus", "M");
305 fitVtxX = histSPDvtxX -> GetFunction("gaus");
306 xMeanVtx = fitVtxX -> GetParameter(1);
307 if (TMath::Abs(xMeanVtx) > 2.) {
309 writeMeanVertexSPD=kTRUE;
312 histSPDvtxY ->Fit("gaus", "M");
313 fitVtxY = histSPDvtxY -> GetFunction("gaus");
314 yMeanVtx = fitVtxY -> GetParameter(1);
315 if (TMath::Abs(yMeanVtx) > 2.) {
317 writeMeanVertexSPD=kTRUE;
320 histSPDvtxZ ->Fit("gaus", "M", "", -12, 12);
321 fitVtxZ = histSPDvtxZ -> GetFunction("gaus");
322 zMeanVtx = fitVtxZ -> GetParameter(1);
323 zSigmaVtx = fitVtxZ -> GetParameter(2);
324 if ((TMath::Abs(zMeanVtx) > 20.) || (zSigmaVtx>12.)) {
325 zMeanVtx = histSPDvtxZ ->GetMean();
326 zSigmaVtx = histSPDvtxZ->GetRMS();
327 writeMeanVertexSPD = kTRUE;
331 else if ((useSPDvtx) && (!spdAvailable)){
332 AliError(Form("Difference between trkVtx and online one, SPD histos not enough entry or SPD 3D vertex off. Writing Mean Vertex SPD"));
333 writeMeanVertexSPD = kTRUE;
337 //check with online position
339 Double_t posOnline[3], sigmaOnline[3];
341 if (useTRKvtx || useITSSAvtx || writeMeanVertexSPD){
342 AliCDBManager *manCheck = AliCDBManager::Instance();
343 manCheck->SetDefaultStorage("raw://");
344 manCheck->SetRun(runNb);
346 AliCDBEntry *entr = manCheck->Get("GRP/Calib/MeanVertexSPD");
348 AliESDVertex *vtxOnline = (AliESDVertex*)entr->GetObject();
350 posOnline[0] = vtxOnline->GetX();
351 posOnline[1] = vtxOnline->GetY();
352 posOnline[2] = vtxOnline->GetZ();
354 sigmaOnline[0] = vtxOnline->GetXRes();
355 sigmaOnline[1] = vtxOnline->GetYRes();
356 sigmaOnline[2] = vtxOnline->GetZRes();
358 //vtxOnline->GetSigmaXYZ(sigmaOnline);
360 if ((TMath::Abs(posOnline[0]-xMeanVtx) > 0.1) || (TMath::Abs(posOnline[1]-yMeanVtx) > 0.1) || (TMath::Abs(posOnline[2]-zMeanVtx) > 1.)){
361 AliWarning(Form("vertex offline far from the online one"));
368 if (writeMeanVertexSPD){
370 AliWarning(Form("Writing Mean Vertex SPD, Mean Vertex not available"));
372 Double_t sigma[3]={0.0150, 0.0150, zSigmaVtx};
374 AliESDVertex *vertex = new AliESDVertex(posOnline, sigma, "vertex");
376 AliCDBId id("GRP/Calib/MeanVertex", runNb, runNb);
378 AliCDBMetaData metaData;
379 metaData.SetBeamPeriod(0); //check!!!!
380 metaData.SetResponsible("Davide Caffarri");
381 metaData.SetComment("Mean Vertex object used in reconstruction");
383 if (!db->Put(vertex, id, &metaData)) {
384 AliError(Form("Error while putting object in storage %s", db->GetURI().Data()));
392 Float_t meanMult = 40.;
394 Float_t resolVtx = 0.05;
396 Double_t xSigmaMult, ySigmaMult, corrXZ, corrYZ, lumiRegSquaredX, lumiRegSquaredY;
397 Double_t covarXZ=0., covarYZ=0.;
399 Bool_t highMultEnvironment = kFALSE;
400 Bool_t highMultppEnvironment = kFALSE;
401 Bool_t lowMultppEnvironment = kFALSE;
403 TF1 *sigmaFitX, *sigmaFitY, *corrFit;
405 //TH1F *histTRKdefMultX=0;
406 //TH1F *histTRKdefMultY=0;
407 TH1F *histTRKHighMultX=0;
408 TH1F *histTRKHighMultY=0;
409 TH2F *histTRKVertexXZ=0;
410 TH2F *histTRKVertexYZ=0;
412 TH2F *histTRKvsMultX=0x0;
413 TH2F *histTRKvsMultY=0x0;
417 //histTRKdefMultX = (TH1F*)list->FindObject("hTRKVertexXdefMult");
418 //histTRKdefMultY = (TH1F*)list->FindObject("hTRKVertexYdefMult");
419 histTRKHighMultX = (TH1F*)list->FindObject("hTRKVertexXHighMult");
420 histTRKHighMultY = (TH1F*)list->FindObject("hTRKVertexYHighMult");
422 histTRKvsMultX = (TH2F*)list->FindObject("hTRKVertexXvsMult");
423 histTRKvsMultY = (TH2F*)list->FindObject("hTRKVertexYvsMult");
425 histTRKVertexXZ = (TH2F*)list->FindObject("hTRKVertexXZ");
426 histTRKVertexYZ = (TH2F*)list->FindObject("hTRKVertexYZ");
430 //histTRKdefMultX = (TH1F*)file->Get("hTRKVertexXdefMult");
431 //histTRKdefMultY = (TH1F*)file->Get("hTRKVertexYdefMult");
432 histTRKHighMultX = (TH1F*)file->Get("hTRKVertexXHighMult");
433 histTRKHighMultY = (TH1F*)file->Get("hTRKVertexYHighMult");
435 histTRKvsMultX = (TH2F*)file->FindObject("hTRKVertexXvsMult");
436 histTRKvsMultY = (TH2F*)file->FindObject("hTRKVertexYvsMult");
438 histTRKVertexXZ = (TH2F*)file->Get("hTRKVertexXZ");
439 histTRKVertexYZ = (TH2F*)file->Get("hTRKVertexYZ");
446 //histTRKdefMultX = (TH1F*)list->FindObject("hITSSAVertexXdefMult");
447 //histTRKdefMultY = (TH1F*)list->FindObject("hITSSAVertexYdefMult");
448 histTRKHighMultX = (TH1F*)list->FindObject("hITSSAVertexXHighMult");
449 histTRKHighMultY = (TH1F*)list->FindObject("hITSSAVertexYHighMult");
451 histTRKvsMultX = (TH2F*)file->FindObject("hITSSAVertexXvsMult");
452 histTRKvsMultY = (TH2F*)file->FindObject("hITSSAVertexYvsMult");
454 histTRKVertexXZ = (TH2F*)list->FindObject("hITSSAVertexXZ");
455 histTRKVertexYZ = (TH2F*)list->FindObject("hITSSAVertexYZ");
459 //histTRKdefMultX = (TH1F*)file->Get("hITSSAVertexXdefMult");
460 //histTRKdefMultY = (TH1F*)file->Get("hITSSAVertexYdefMult");
461 histTRKHighMultX = (TH1F*)file->Get("hITSSAVertexXHighMult");
462 histTRKHighMultY = (TH1F*)file->Get("hITSSAVertexYHighMult");
464 histTRKvsMultX = (TH2F*)file->FindObject("hITSSAVertexXvsMult");
465 histTRKvsMultY = (TH2F*)file->FindObject("hITSSAVertexYvsMult");
467 histTRKVertexXZ = (TH2F*)file->Get("hITSSAVertexXZ");
468 histTRKVertexYZ = (TH2F*)file->Get("hITSSAVertexYZ");
476 Float_t nEntriesMultX=0, nEntriesMultY=0.;
478 if ((histTRKHighMultX) && (histTRKHighMultY)){
480 projXvsMult = (TH1D*)histTRKvsMultX->ProjectionX("projXHighMultPbPb", 150, 300);
481 projYvsMult = (TH1D*)histTRKvsMultY->ProjectionX("projYHighMultPbPb", 150, 300);
483 nEntriesMultX = projXvsMult->GetEffectiveEntries();
484 nEntriesMultY = projYvsMult->GetEffectiveEntries();
486 if ((nEntriesMultX >100) && (nEntriesMultY>100)) {
487 AliWarning(Form("Setting High Mulitplicity environment"));
488 highMultEnvironment = kTRUE;
492 projXvsMult = (TH1D*)histTRKvsMultX->ProjectionX("projXHighMultPbPb", 10, 30);
493 projYvsMult = (TH1D*)histTRKvsMultY->ProjectionX("projYHighMultPbPb", 10, 30);
495 nEntriesMultX = projXvsMult->GetEffectiveEntries();
496 nEntriesMultY = projYvsMult->GetEffectiveEntries();
498 if ((nEntriesMultX >100) && (nEntriesMultY>100)) {
499 AliWarning(Form("Setting high pp Mulitplicity environment or p-A high multiplicity"));
500 highMultppEnvironment=kTRUE;
504 projXvsMult = (TH1D*)histTRKvsMultX->ProjectionX("projXHighMultPbPb", 3, 5);
505 projYvsMult = (TH1D*)histTRKvsMultY->ProjectionX("projYHighMultPbPb", 3, 5);
507 nEntriesMultX = projXvsMult->GetEffectiveEntries();
508 nEntriesMultY = projYvsMult->GetEffectiveEntries();
510 if ((nEntriesMultX >100) && (nEntriesMultY>100)) {
511 AliWarning(Form("Setting low pp Mulitplicity environment"));
512 lowMultppEnvironment=kTRUE;
520 if (lowMultppEnvironment==kTRUE) {
522 if ((projXvsMult->GetEntries() < 40.) || (projYvsMult->GetEntries() < 40.)){
523 AliWarning(Form("histos for lumi reg calculation not found, default value set"));
526 fStatus=kLumiRegCovMatrixProblem;
529 projXvsMult -> Fit("gaus", "M", "", -0.4, 0.4);
530 sigmaFitX = projXvsMult -> GetFunction("gaus");
531 xSigmaMult = sigmaFitX->GetParameter(2);
533 lumiRegSquaredX = (xSigmaMult*xSigmaMult - ((resolVtx*resolVtx)/TMath::Power(meanMult, p2)));
534 if (lumiRegSquaredX < 0 || lumiRegSquaredX < 1E-5) {
535 AliWarning(Form("Difficult luminous region determination X, keep convoluted sigma"));
536 xSigmaVtx = xSigmaMult;
537 fStatus=kLumiRegCovMatrixProblem;
541 if (lumiRegSquaredX > 0 && lumiRegSquaredX < 0.0005){
542 xSigmaVtx = TMath::Sqrt(lumiRegSquaredX);
543 xSigmaVtx = xSigmaVtx*1.1;
547 AliWarning(Form("Not possible to define a luminous region X. Default values set"));
549 fStatus=kLumiRegCovMatrixProblem;
553 projYvsMult -> Fit("gaus", "M", "", -0.2, 0.6);
554 sigmaFitY = projYvsMult -> GetFunction("gaus");
555 ySigmaMult = sigmaFitY->GetParameter(2);
557 lumiRegSquaredY= (ySigmaMult*ySigmaMult - ((resolVtx*resolVtx)/TMath::Power(meanMult, p2)));
559 if (lumiRegSquaredY < 0 || lumiRegSquaredY < 1E-5) {
560 AliWarning(Form("Difficult luminous region determination Y, keep convoluted sigma"));
561 ySigmaVtx = ySigmaMult;
562 fStatus=kLumiRegCovMatrixProblem;
566 if (lumiRegSquaredY > 0 && lumiRegSquaredY < 0.0005){
567 ySigmaVtx = TMath::Sqrt(lumiRegSquaredY);
568 ySigmaVtx = ySigmaVtx*1.1;
571 AliWarning(Form("Not possible to define a luminous region Y. Default values set"));
573 fStatus=kLumiRegCovMatrixProblem;
578 TProfile *htrkXZ = histTRKVertexXZ ->ProfileY();
579 htrkXZ -> Fit("pol1", "M", "", -10., 10.);
580 corrFit = htrkXZ->GetFunction("pol1");
581 corrXZ = corrFit->GetParameter(1);
583 if (TMath::Abs(corrXZ) > 0.01) {
584 AliWarning(Form("Problems in the correlation fitting, not update the covariance matrix"));
586 fStatus=kLumiRegCovMatrixProblem;
589 covarXZ = corrXZ * zSigmaVtx*zSigmaVtx;
593 TProfile *htrkYZ = histTRKVertexYZ ->ProfileY();
594 htrkYZ -> Fit("pol1", "M", "", -10., 10.);
595 corrFit = htrkYZ->GetFunction("pol1");
596 corrYZ = corrFit->GetParameter(1);
598 if (TMath::Abs(corrYZ) > 0.01) {
599 AliWarning(Form("Problems in the correlation fitting, not update the covariance matrix"));
601 fStatus=kLumiRegCovMatrixProblem;
604 covarYZ = corrYZ*zSigmaVtx*zSigmaVtx;
610 if (highMultEnvironment==kTRUE || highMultppEnvironment==kTRUE){
612 projXvsMult -> Fit("gaus", "M", "", -0.4, 0.4);
613 sigmaFitX = projXvsMult -> GetFunction("gaus");
614 xSigmaMult = sigmaFitX->GetParameter(2);
616 if ((xSigmaMult <0) || (xSigmaMult>0.03)){
617 AliWarning(Form("Problems with luminosiy region determination, update of the postion only"));
620 fStatus=kLumiRegCovMatrixProblem;
623 xSigmaVtx = xSigmaMult;
624 xSigmaVtx = xSigmaVtx*1.1;
627 projYvsMult -> Fit("gaus", "M", "", -0.2, 0.5);
628 TCanvas *c = new TCanvas("nwC", "nwC");
631 sigmaFitY = projYvsMult -> GetFunction("gaus");
632 ySigmaMult = sigmaFitY->GetParameter(2);
634 if ((ySigmaMult <0) || (ySigmaMult>0.03)){
635 AliWarning(Form("Problems with luminosiy region determination, update of the postion only"));
638 fStatus=kLumiRegCovMatrixProblem;
641 ySigmaVtx = ySigmaMult;
642 ySigmaVtx = ySigmaVtx*1.1;
645 TProfile *htrkXZ = histTRKVertexXZ ->ProfileY();
646 htrkXZ -> Fit("pol1", "M", "", -10., 10.);
647 corrFit = htrkXZ->GetFunction("pol1");
648 corrXZ = corrFit->GetParameter(1);
650 if (TMath::Abs(corrXZ) > 0.01) {
651 AliWarning(Form("Problems in the correlation fitting, not update the covariance matrix"));
653 fStatus=kLumiRegCovMatrixProblem;
656 covarXZ = corrXZ * zSigmaVtx*zSigmaVtx;
660 TProfile *htrkYZ = histTRKVertexYZ ->ProfileY();
661 htrkYZ -> Fit("pol1", "M", "", -10., 10.);
662 corrFit = htrkYZ->GetFunction("pol1");
663 corrYZ = corrFit->GetParameter(1);
665 if (TMath::Abs(corrYZ) > 0.01) {
666 AliWarning(Form("Problems in the correlation fitting, not update the covariance matrix"));
668 fStatus=kLumiRegCovMatrixProblem;
671 covarYZ = corrYZ*zSigmaVtx*zSigmaVtx;
676 Double_t position[3], covMatrix[6];
680 position[0] = xMeanVtx;
681 position[1] = yMeanVtx;
682 position[2] = zMeanVtx;
684 covMatrix[0] = xSigmaVtx*xSigmaVtx;
685 covMatrix[1] = 0.; //xy
686 covMatrix[2] = ySigmaVtx*ySigmaVtx;
687 covMatrix[3] = covarXZ;
688 covMatrix[4] = covarYZ;
689 covMatrix[5] = zSigmaVtx*zSigmaVtx;
692 //Printf ("%f, %f, %f, %f", xSigmaVtx, ySigmaVtx, covarXZ, covarYZ);
694 AliESDVertex *vertex = new AliESDVertex(position, covMatrix, chi2, nContr, "vertex");
696 AliCDBId id("GRP/Calib/MeanVertex", runNb, runNb);
698 AliCDBMetaData metaData;
699 metaData.SetBeamPeriod(0); //check!!!!
700 metaData.SetResponsible("Davide Caffarri");
701 metaData.SetComment("Mean Vertex object used in reconstruction");
703 if (!db->Put(vertex, id, &metaData)) {
704 AliError(Form("Error while putting object in storage %s", db->GetURI().Data()));
710 Int_t status=GetStatus();
712 AliInfo(Form("MeanVertex calibration successful: %s (status=%d)", fgkStatusCodeName[fStatus], status));
714 else if (status > 0) {
715 AliInfo(Form("MeanVertex calibration failed: %s (status=%d)", fgkStatusCodeName[fStatus], status));
717 else if (status < 0) {
718 AliInfo(Form("MeanVertex calibration but not fatal error: %s (status=%d)", fgkStatusCodeName[fStatus], status));
724 //__________________________________________________________________________
725 Int_t AliMeanVertexPreprocessorOffline::GetStatus(){
732 /* OK, return zero */
737 /* non-fatal error, return negative status */
739 case kWriteMeanVertexSPD:
740 case kUseOfflineSPDvtx:
741 case kLumiRegCovMatrixProblem:
745 /* fatal error, return positive status */
751 /* anything else, return negative large number */
757 /* should never arrive here, anyway return negative large number */