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;
299 AliDebug(2, Form("xHistoRMS = %f, yHistoRMS = %f, zHistoRMS = %f", xHistoRMS, yHistoRMS, zHistoRMS));
303 if ((useSPDvtx) && (spdAvailable) && (!vertexerSPD3Doff)){
305 histSPDvtxX ->Fit("gaus", "M");
306 fitVtxX = histSPDvtxX -> GetFunction("gaus");
307 xMeanVtx = fitVtxX -> GetParameter(1);
308 xSigmaVtx = fitVtxX -> GetParameter(2);
309 if (TMath::Abs(xMeanVtx) > 2.) {
311 writeMeanVertexSPD=kTRUE;
314 histSPDvtxY ->Fit("gaus", "M");
315 fitVtxY = histSPDvtxY -> GetFunction("gaus");
316 yMeanVtx = fitVtxY -> GetParameter(1);
317 ySigmaVtx = fitVtxY -> GetParameter(2);
318 if (TMath::Abs(yMeanVtx) > 2.) {
320 writeMeanVertexSPD=kTRUE;
323 histSPDvtxZ ->Fit("gaus", "M", "", -12, 12);
324 fitVtxZ = histSPDvtxZ -> GetFunction("gaus");
325 zMeanVtx = fitVtxZ -> GetParameter(1);
326 zSigmaVtx = fitVtxZ -> GetParameter(2);
327 if ((TMath::Abs(zMeanVtx) > 20.) || (zSigmaVtx>12.)) {
328 zMeanVtx = histSPDvtxZ ->GetMean();
329 zSigmaVtx = histSPDvtxZ->GetRMS();
330 writeMeanVertexSPD = kTRUE;
334 else if ((useSPDvtx) && (!spdAvailable)){
335 AliError(Form("Difference between trkVtx and online one, SPD histos not enough entry or SPD 3D vertex off. Writing Mean Vertex SPD"));
336 writeMeanVertexSPD = kTRUE;
340 //check with online position
342 Double_t posOnline[3], sigmaOnline[3];
344 if (useTRKvtx || useITSSAvtx || writeMeanVertexSPD){
345 AliCDBManager *manCheck = AliCDBManager::Instance();
346 manCheck->SetDefaultStorage("raw://");
347 manCheck->SetRun(runNb);
349 AliCDBEntry *entr = manCheck->Get("GRP/Calib/MeanVertexSPD");
351 AliESDVertex *vtxOnline = (AliESDVertex*)entr->GetObject();
353 posOnline[0] = vtxOnline->GetX();
354 posOnline[1] = vtxOnline->GetY();
355 posOnline[2] = vtxOnline->GetZ();
357 sigmaOnline[0] = vtxOnline->GetXRes();
358 sigmaOnline[1] = vtxOnline->GetYRes();
359 sigmaOnline[2] = vtxOnline->GetZRes();
361 AliDebug(2, Form("sigmaOnline[0] = %f, sigmaOnline[1] = %f, sigmaOnline[2] = %f", sigmaOnline[0], sigmaOnline[1], sigmaOnline[2]));
362 //vtxOnline->GetSigmaXYZ(sigmaOnline);
364 if ((TMath::Abs(posOnline[0]-xMeanVtx) > 0.1) || (TMath::Abs(posOnline[1]-yMeanVtx) > 0.1) || (TMath::Abs(posOnline[2]-zMeanVtx) > 1.)){
365 AliWarning(Form("vertex offline far from the online one"));
372 if (writeMeanVertexSPD){
374 AliWarning(Form("Writing Mean Vertex SPD, Mean Vertex not available"));
376 Double_t sigma[3]={0.0150, 0.0150, zSigmaVtx};
378 AliESDVertex *vertex = new AliESDVertex(posOnline, sigma, "vertex");
380 AliCDBId id("GRP/Calib/MeanVertex", runNb, runNb);
382 AliCDBMetaData metaData;
383 metaData.SetBeamPeriod(0); //check!!!!
384 metaData.SetResponsible("Davide Caffarri");
385 metaData.SetComment("Mean Vertex object used in reconstruction");
387 if (!db->Put(vertex, id, &metaData)) {
388 AliError(Form("Error while putting object in storage %s", db->GetURI().Data()));
396 Float_t meanMult = 40.;
398 Float_t resolVtx = 0.05;
400 Double_t xSigmaMult, ySigmaMult, corrXZ, corrYZ, lumiRegSquaredX, lumiRegSquaredY;
401 Double_t covarXZ=0., covarYZ=0.;
403 Bool_t highMultEnvironment = kFALSE;
404 Bool_t highMultppEnvironment = kFALSE;
405 Bool_t lowMultppEnvironment = kFALSE;
407 TF1 *sigmaFitX, *sigmaFitY, *corrFit;
409 //TH1F *histTRKdefMultX=0;
410 //TH1F *histTRKdefMultY=0;
411 TH1F *histTRKHighMultX=0;
412 TH1F *histTRKHighMultY=0;
413 TH2F *histTRKVertexXZ=0;
414 TH2F *histTRKVertexYZ=0;
416 TH2F *histTRKvsMultX=0x0;
417 TH2F *histTRKvsMultY=0x0;
421 //histTRKdefMultX = (TH1F*)list->FindObject("hTRKVertexXdefMult");
422 //histTRKdefMultY = (TH1F*)list->FindObject("hTRKVertexYdefMult");
423 histTRKHighMultX = (TH1F*)list->FindObject("hTRKVertexXHighMult");
424 histTRKHighMultY = (TH1F*)list->FindObject("hTRKVertexYHighMult");
426 histTRKvsMultX = (TH2F*)list->FindObject("hTRKVertexXvsMult");
427 histTRKvsMultY = (TH2F*)list->FindObject("hTRKVertexYvsMult");
429 histTRKVertexXZ = (TH2F*)list->FindObject("hTRKVertexXZ");
430 histTRKVertexYZ = (TH2F*)list->FindObject("hTRKVertexYZ");
434 //histTRKdefMultX = (TH1F*)file->Get("hTRKVertexXdefMult");
435 //histTRKdefMultY = (TH1F*)file->Get("hTRKVertexYdefMult");
436 histTRKHighMultX = (TH1F*)file->Get("hTRKVertexXHighMult");
437 histTRKHighMultY = (TH1F*)file->Get("hTRKVertexYHighMult");
439 histTRKvsMultX = (TH2F*)file->FindObject("hTRKVertexXvsMult");
440 histTRKvsMultY = (TH2F*)file->FindObject("hTRKVertexYvsMult");
442 histTRKVertexXZ = (TH2F*)file->Get("hTRKVertexXZ");
443 histTRKVertexYZ = (TH2F*)file->Get("hTRKVertexYZ");
450 //histTRKdefMultX = (TH1F*)list->FindObject("hITSSAVertexXdefMult");
451 //histTRKdefMultY = (TH1F*)list->FindObject("hITSSAVertexYdefMult");
452 histTRKHighMultX = (TH1F*)list->FindObject("hITSSAVertexXHighMult");
453 histTRKHighMultY = (TH1F*)list->FindObject("hITSSAVertexYHighMult");
455 histTRKvsMultX = (TH2F*)file->FindObject("hITSSAVertexXvsMult");
456 histTRKvsMultY = (TH2F*)file->FindObject("hITSSAVertexYvsMult");
458 histTRKVertexXZ = (TH2F*)list->FindObject("hITSSAVertexXZ");
459 histTRKVertexYZ = (TH2F*)list->FindObject("hITSSAVertexYZ");
463 //histTRKdefMultX = (TH1F*)file->Get("hITSSAVertexXdefMult");
464 //histTRKdefMultY = (TH1F*)file->Get("hITSSAVertexYdefMult");
465 histTRKHighMultX = (TH1F*)file->Get("hITSSAVertexXHighMult");
466 histTRKHighMultY = (TH1F*)file->Get("hITSSAVertexYHighMult");
468 histTRKvsMultX = (TH2F*)file->FindObject("hITSSAVertexXvsMult");
469 histTRKvsMultY = (TH2F*)file->FindObject("hITSSAVertexYvsMult");
471 histTRKVertexXZ = (TH2F*)file->Get("hITSSAVertexXZ");
472 histTRKVertexYZ = (TH2F*)file->Get("hITSSAVertexYZ");
480 Float_t nEntriesMultX=0, nEntriesMultY=0.;
482 if ((histTRKHighMultX) && (histTRKHighMultY)){
484 projXvsMult = (TH1D*)histTRKvsMultX->ProjectionX("projXHighMultPbPb", 150, 300);
485 projYvsMult = (TH1D*)histTRKvsMultY->ProjectionX("projYHighMultPbPb", 150, 300);
487 nEntriesMultX = projXvsMult->GetEffectiveEntries();
488 nEntriesMultY = projYvsMult->GetEffectiveEntries();
490 if ((nEntriesMultX >100) && (nEntriesMultY>100)) {
491 AliWarning(Form("Setting High Mulitplicity environment"));
492 highMultEnvironment = kTRUE;
496 projXvsMult = (TH1D*)histTRKvsMultX->ProjectionX("projXHighMultPbPb", 10, 30);
497 projYvsMult = (TH1D*)histTRKvsMultY->ProjectionX("projYHighMultPbPb", 10, 30);
499 nEntriesMultX = projXvsMult->GetEffectiveEntries();
500 nEntriesMultY = projYvsMult->GetEffectiveEntries();
502 if ((nEntriesMultX >100) && (nEntriesMultY>100)) {
503 AliWarning(Form("Setting high pp Mulitplicity environment or p-A high multiplicity"));
504 highMultppEnvironment=kTRUE;
508 projXvsMult = (TH1D*)histTRKvsMultX->ProjectionX("projXHighMultPbPb", 3, 5);
509 projYvsMult = (TH1D*)histTRKvsMultY->ProjectionX("projYHighMultPbPb", 3, 5);
511 nEntriesMultX = projXvsMult->GetEffectiveEntries();
512 nEntriesMultY = projYvsMult->GetEffectiveEntries();
514 if ((nEntriesMultX >100) && (nEntriesMultY>100)) {
515 AliWarning(Form("Setting low pp Mulitplicity environment"));
516 lowMultppEnvironment=kTRUE;
524 if (lowMultppEnvironment==kTRUE) {
526 if ((projXvsMult->GetEntries() < 40.) || (projYvsMult->GetEntries() < 40.)){
527 AliWarning(Form("histos for lumi reg calculation not found, default value set"));
530 fStatus=kLumiRegCovMatrixProblem;
533 projXvsMult -> Fit("gaus", "M", "", -0.4, 0.4);
534 sigmaFitX = projXvsMult -> GetFunction("gaus");
535 xSigmaMult = sigmaFitX->GetParameter(2);
537 lumiRegSquaredX = (xSigmaMult*xSigmaMult - ((resolVtx*resolVtx)/TMath::Power(meanMult, p2)));
538 if (lumiRegSquaredX < 0 || lumiRegSquaredX < 1E-5) {
539 AliWarning(Form("Difficult luminous region determination X, keep convoluted sigma"));
540 xSigmaVtx = xSigmaMult;
541 fStatus=kLumiRegCovMatrixProblem;
545 if (lumiRegSquaredX > 0 && lumiRegSquaredX < 0.0005){
546 xSigmaVtx = TMath::Sqrt(lumiRegSquaredX);
547 xSigmaVtx = xSigmaVtx*1.1;
551 AliWarning(Form("Not possible to define a luminous region X. Default values set"));
553 fStatus=kLumiRegCovMatrixProblem;
557 projYvsMult -> Fit("gaus", "M", "", -0.2, 0.6);
558 sigmaFitY = projYvsMult -> GetFunction("gaus");
559 ySigmaMult = sigmaFitY->GetParameter(2);
561 lumiRegSquaredY= (ySigmaMult*ySigmaMult - ((resolVtx*resolVtx)/TMath::Power(meanMult, p2)));
563 if (lumiRegSquaredY < 0 || lumiRegSquaredY < 1E-5) {
564 AliWarning(Form("Difficult luminous region determination Y, keep convoluted sigma"));
565 ySigmaVtx = ySigmaMult;
566 fStatus=kLumiRegCovMatrixProblem;
570 if (lumiRegSquaredY > 0 && lumiRegSquaredY < 0.0005){
571 ySigmaVtx = TMath::Sqrt(lumiRegSquaredY);
572 ySigmaVtx = ySigmaVtx*1.1;
575 AliWarning(Form("Not possible to define a luminous region Y. Default values set"));
577 fStatus=kLumiRegCovMatrixProblem;
582 TProfile *htrkXZ = histTRKVertexXZ ->ProfileY();
583 htrkXZ -> Fit("pol1", "M", "", -10., 10.);
584 corrFit = htrkXZ->GetFunction("pol1");
585 corrXZ = corrFit->GetParameter(1);
587 if (TMath::Abs(corrXZ) > 0.01) {
588 AliWarning(Form("Problems in the correlation fitting, not update the covariance matrix"));
590 fStatus=kLumiRegCovMatrixProblem;
593 covarXZ = corrXZ * zSigmaVtx*zSigmaVtx;
597 TProfile *htrkYZ = histTRKVertexYZ ->ProfileY();
598 htrkYZ -> Fit("pol1", "M", "", -10., 10.);
599 corrFit = htrkYZ->GetFunction("pol1");
600 corrYZ = corrFit->GetParameter(1);
602 if (TMath::Abs(corrYZ) > 0.01) {
603 AliWarning(Form("Problems in the correlation fitting, not update the covariance matrix"));
605 fStatus=kLumiRegCovMatrixProblem;
608 covarYZ = corrYZ*zSigmaVtx*zSigmaVtx;
614 if (highMultEnvironment==kTRUE || highMultppEnvironment==kTRUE){
616 projXvsMult -> Fit("gaus", "M", "", -0.4, 0.4);
617 sigmaFitX = projXvsMult -> GetFunction("gaus");
618 xSigmaMult = sigmaFitX->GetParameter(2);
620 if ((xSigmaMult <0) || (xSigmaMult>0.03)){
621 AliWarning(Form("Problems with luminosiy region determination, update of the postion only"));
624 fStatus=kLumiRegCovMatrixProblem;
627 xSigmaVtx = xSigmaMult;
628 xSigmaVtx = xSigmaVtx*1.1;
631 projYvsMult -> Fit("gaus", "M", "", -0.2, 0.5);
632 TCanvas *c = new TCanvas("nwC", "nwC");
635 sigmaFitY = projYvsMult -> GetFunction("gaus");
636 ySigmaMult = sigmaFitY->GetParameter(2);
638 if ((ySigmaMult <0) || (ySigmaMult>0.03)){
639 AliWarning(Form("Problems with luminosiy region determination, update of the postion only"));
642 fStatus=kLumiRegCovMatrixProblem;
645 ySigmaVtx = ySigmaMult;
646 ySigmaVtx = ySigmaVtx*1.1;
649 TProfile *htrkXZ = histTRKVertexXZ ->ProfileY();
650 htrkXZ -> Fit("pol1", "M", "", -10., 10.);
651 corrFit = htrkXZ->GetFunction("pol1");
652 corrXZ = corrFit->GetParameter(1);
654 if (TMath::Abs(corrXZ) > 0.01) {
655 AliWarning(Form("Problems in the correlation fitting, not update the covariance matrix"));
657 fStatus=kLumiRegCovMatrixProblem;
660 covarXZ = corrXZ * zSigmaVtx*zSigmaVtx;
664 TProfile *htrkYZ = histTRKVertexYZ ->ProfileY();
665 htrkYZ -> Fit("pol1", "M", "", -10., 10.);
666 corrFit = htrkYZ->GetFunction("pol1");
667 corrYZ = corrFit->GetParameter(1);
669 if (TMath::Abs(corrYZ) > 0.01) {
670 AliWarning(Form("Problems in the correlation fitting, not update the covariance matrix"));
672 fStatus=kLumiRegCovMatrixProblem;
675 covarYZ = corrYZ*zSigmaVtx*zSigmaVtx;
680 Double_t position[3], covMatrix[6];
684 position[0] = xMeanVtx;
685 position[1] = yMeanVtx;
686 position[2] = zMeanVtx;
688 covMatrix[0] = xSigmaVtx*xSigmaVtx;
689 covMatrix[1] = 0.; //xy
690 covMatrix[2] = ySigmaVtx*ySigmaVtx;
691 covMatrix[3] = covarXZ;
692 covMatrix[4] = covarYZ;
693 covMatrix[5] = zSigmaVtx*zSigmaVtx;
696 //Printf ("%f, %f, %f, %f", xSigmaVtx, ySigmaVtx, covarXZ, covarYZ);
698 AliESDVertex *vertex = new AliESDVertex(position, covMatrix, chi2, nContr, "vertex");
700 AliCDBId id("GRP/Calib/MeanVertex", runNb, runNb);
702 AliCDBMetaData metaData;
703 metaData.SetBeamPeriod(0); //check!!!!
704 metaData.SetResponsible("Davide Caffarri");
705 metaData.SetComment("Mean Vertex object used in reconstruction");
707 if (!db->Put(vertex, id, &metaData)) {
708 AliError(Form("Error while putting object in storage %s", db->GetURI().Data()));
714 Int_t status=GetStatus();
716 AliInfo(Form("MeanVertex calibration successful: %s (status=%d)", fgkStatusCodeName[fStatus], status));
718 else if (status > 0) {
719 AliInfo(Form("MeanVertex calibration failed: %s (status=%d)", fgkStatusCodeName[fStatus], status));
721 else if (status < 0) {
722 AliInfo(Form("MeanVertex calibration but not fatal error: %s (status=%d)", fgkStatusCodeName[fStatus], status));
728 //__________________________________________________________________________
729 Int_t AliMeanVertexPreprocessorOffline::GetStatus(){
736 /* OK, return zero */
741 /* non-fatal error, return negative status */
743 case kWriteMeanVertexSPD:
744 case kUseOfflineSPDvtx:
745 case kLumiRegCovMatrixProblem:
749 /* fatal error, return positive status */
755 /* anything else, return negative large number */
761 /* should never arrive here, anyway return negative large number */