/************************************************************************** * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * * * * Author: The ALICE Off-line Project. * * Contributors are mentioned in the code where appropriate. * * * * Permission to use, copy, modify and distribute this software and its * * documentation strictly for non-commercial purposes is hereby granted * * without fee, provided that the above copyright notice appears in all * * copies and that both the copyright notice and this permission notice * * appear in the supporting documentation. The authors make no claims * * about the suitability of this software for any purpose. It is * * provided "as is" without express or implied warranty. * **************************************************************************/ //////////////////////////////////////////////////////////////////////////// // // // Gain calibration using tracks // // The main goal: // 1.) Inner TPC gain alignement - (parabolic) parameterization (inside of one sector) // 2.) Angular and z-position correction (parabolic) parameterization // 3.) Sector gain alignment // // Following histograms are accumulated // a.) Simple 1D histograms per chamber // b.) Profile histograms per chamber - local x dependence // c.) 2D Profile histograms - local x - fi dependence // // To get the gain map - the simple solution - use the histograms - is not enough // The resulting mean amplitude map depends strongly on the track topology // These dependence can be reduced, taking into account angular effect, and diffusion effect // Using proper fit modeles // // // // // === Calibration class for gain calibration using tracks === // // 1) Genereal idea // ================ // A 6-parametric parabolic function // // G(x, y) = p0 + p1*x + p2*y + p3*x^2 + p4*y^2 + p5 * x*y // // is fitted to the maximum charge values or total charge values of // all the clusters contained in the tracks that are added to this // object. This fit is performed for each read out chamber, in fact even // for each type of pad sizes (thus for one segment, which consists of // an IROC and an OROC, there are three fitters used, corresponding to // the three pad sizes). The coordinate origin is at the center of the // particular pad size region on each ROC. // // Because of the Landau nature of the charge deposition we use // different "types" of fitters instead of one to minimize the effect // of the long Landau tail. The difference between the fitters is only // the charge value, that is put into them, i.e. the charge is subject // to a transformation. At this point we use three different fit types: // // a) simple: the charge is put in as it is // b) sqrt: the square root of the charge is put into the fitter // c) log: fgkM * Log(1+q/fgkM) is put into the fitter, with // q being the untransformed charge and fgkM=25 // // The results of the fits may be visualized and further used by // creating an AliTPCCalROC or AliTPCCalPad. You may specify to undo // the transformation and/or to normalize to the pad size. // // Not every track you add to this object is actually used for // calibration. There are some cuts and conditions to exclude bad // tracks, e.g. a pt cut to cut out tracks with too much charge // deposition or a cut on edge clusters which are not fully // registered and don't give a usable signal. // // 2) Interface / usage // ==================== // For each track to be added you need to call Process(). // This method expects an AliTPCseed, which contains the necessary // cluster information. At the moment of writing this information // is stored in an AliESDfriend corresponding to an AliESD. // You may also call AddTrack() if you don't want the cuts and // other quality conditions to kick in (thus forcing the object to // accept the track) or AddCluster() for adding single clusters. // Call one of the Evaluate functions to evaluate the fitter(s) and // to retrieve the fit parameters, erros and so on. You can also // do this later on by using the different Getters. // // The visualization methods CreateFitCalPad() and CreateFitCalROC() // are straight forward to use. // // Note: If you plan to write this object to a ROOT file, make sure // you evaluate all the fitters *before* writing, because due // to a bug in the fitter component writing fitters doesn't // work properly (yet). Be aware that you cannot re-evaluate // the fitters after loading this object from file. // (This will be gone for a new ROOT version > v5-17-05) // //////////////////////////////////////////////////////////////////////////// #include #include #include "TSystem.h" #include "TMatrixD.h" #include "TTreeStream.h" #include "TF1.h" #include "AliTPCParamSR.h" #include "AliTPCClusterParam.h" #include "AliTrackPointArray.h" #include "TCint.h" #include "AliTPCcalibTracksGain.h" #include #include #include #include #include #include #include // // AliRoot includes // #include "AliMagF.h" #include "AliMathBase.h" // #include "AliTPCROC.h" #include "AliTPCParamSR.h" #include "AliTPCCalROC.h" #include "AliTPCCalPad.h" #include "AliTPCClusterParam.h" // #include "AliTracker.h" #include "AliESD.h" #include "AliESDtrack.h" #include "AliESDfriend.h" #include "AliESDfriendTrack.h" #include "AliTPCseed.h" #include "AliTPCclusterMI.h" #include "AliTPCcalibTracksCuts.h" #include "AliTPCFitPad.h" // REMOVE ALL OF THIS #include #include "AliESDEvent.h" /* TFile f("TPCCalibTracksGain.root") gSystem->Load("libPWG1.so") AliTreeDraw comp comp.SetTree(dEdx) Double_t chi2; TVectorD vec(3) TMatrixD mat(3,3) TString * str = comp.FitPlane("Cl.fQ/dedxQ.fElements[0]","Cl.fY++Cl.fX","Cl.fDetector<36",chi2,vec,mat) */ ClassImp(AliTPCcalibTracksGain) const Bool_t AliTPCcalibTracksGain::fgkUseTotalCharge = kTRUE; const Double_t AliTPCcalibTracksGain::fgkM = 25.; const char* AliTPCcalibTracksGain::fgkDebugStreamFileName = "TPCCalibTracksGain.root"; AliTPCParamSR* AliTPCcalibTracksGain::fgTPCparam = new AliTPCParamSR(); AliTPCcalibTracksGain::AliTPCcalibTracksGain() : AliTPCcalibBase(), fDebugStream(0), //! debug stream for debugging fCuts(0), // cuts that are used for sieving the tracks used for calibration // // Simple Array of histograms // fArrayQM(0), // Qmax normalized fArrayQT(0), // Qtot normalized fProfileArrayQM(0), // Qmax normalized versus local X fProfileArrayQT(0), // Qtot normalized versus local X fProfileArrayQM2D(0), // Qmax normalized versus local X and phi fProfileArrayQT2D(0), // Qtot normalized versus local X and phi // // Fitters // fSimpleFitter(0), // simple fitter for short pads fSqrtFitter(0), // sqrt fitter for medium pads fLogFitter(0), // log fitter for long pads fFitter0M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter1M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter2M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter0T(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter1T(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter2T(0), // fitting of the atenuation, angular correction, and mean chamber gain fSingleSectorFitter(0), // just for debugging // // Counters // fTotalTracks(0), // just for debugging fAcceptedTracks(0), // just for debugging fDebugCalPadRaw(0), // just for debugging fDebugCalPadCorr(0), // just for debugging fPrevIter(0) // the calibration object in its previous iteration (will not be owned by the new object, don't forget to delete it!) { // // Default constructor. // } AliTPCcalibTracksGain::AliTPCcalibTracksGain(const AliTPCcalibTracksGain& obj) : AliTPCcalibBase(obj), fDebugStream(0), //! debug stream for debugging fCuts(obj.fCuts), // cuts that are used for sieving the tracks used for calibration fArrayQM(0), // Qmax normalized fArrayQT(0), // Qtot normalized // // Simple Histograms // fProfileArrayQM(obj.fProfileArrayQM), // Qmax normalized versus local X fProfileArrayQT(obj.fProfileArrayQT), // Qtot normalized versus local X fProfileArrayQM2D(obj.fProfileArrayQM2D), // Qmax normalized versus local X and phi fProfileArrayQT2D(obj.fProfileArrayQT2D), // Qtot normalized versus local X and phi // // Fitters // fSimpleFitter(obj.fSimpleFitter), // simple fitter for short pads fSqrtFitter(obj.fSqrtFitter), // sqrt fitter for medium pads fLogFitter(obj.fLogFitter), // log fitter for long pads fFitter0M(obj.fFitter0M), fFitter1M(obj.fFitter1M), fFitter2M(obj.fFitter2M), fFitter0T(obj.fFitter0T), fFitter1T(obj.fFitter1T), fFitter2T(obj.fFitter2T), fSingleSectorFitter(obj.fSingleSectorFitter), // just for debugging // // Counters // fTotalTracks(obj.fTotalTracks), // just for debugging fAcceptedTracks(obj.fAcceptedTracks), // just for debugging fDebugCalPadRaw(obj.fDebugCalPadRaw), // just for debugging fDebugCalPadCorr(obj.fDebugCalPadCorr), // just for debugging fPrevIter(obj.fPrevIter) // the calibration object in its previous iteration (will not be owned by the new object, don't forget to delete it!) { // // Copy constructor. // } AliTPCcalibTracksGain& AliTPCcalibTracksGain::operator=(const AliTPCcalibTracksGain& rhs) { // // Assignment operator. // if (this != &rhs) { TNamed::operator=(rhs); fDebugCalPadRaw = new AliTPCCalPad(*(rhs.fDebugCalPadRaw)); fDebugCalPadCorr = new AliTPCCalPad(*(rhs.fDebugCalPadCorr)); fSimpleFitter = new AliTPCFitPad(*(rhs.fSimpleFitter)); fSqrtFitter = new AliTPCFitPad(*(rhs.fSqrtFitter)); fLogFitter = new AliTPCFitPad(*(rhs.fLogFitter)); fSingleSectorFitter = new AliTPCFitPad(*(rhs.fSingleSectorFitter)); fPrevIter = new AliTPCcalibTracksGain(*(rhs.fPrevIter)); fCuts = new AliTPCcalibTracksCuts(*(rhs.fCuts)); } return *this; } AliTPCcalibTracksGain::AliTPCcalibTracksGain(const char* name, const char* title, AliTPCcalibTracksCuts* cuts, TNamed* /*debugStreamPrefix*/, AliTPCcalibTracksGain* prevIter) : AliTPCcalibBase(), fDebugStream(0), //! debug stream for debugging fCuts(0), // cuts that are used for sieving the tracks used for calibration fArrayQM(0), // Qmax normalized fArrayQT(0), // Qtot normalized // // Simple Histograms // fProfileArrayQM(0), // Qmax normalized versus local X fProfileArrayQT(0), // Qtot normalized versus local X fProfileArrayQM2D(0), // Qmax normalized versus local X and phi fProfileArrayQT2D(0), // Qtot normalized versus local X and phi // // Fitters // fSimpleFitter(0), // simple fitter for short pads fSqrtFitter(0), // sqrt fitter for medium pads fLogFitter(0), // log fitter for long pads fFitter0M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter1M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter2M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter0T(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter1T(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter2T(0), // fitting of the atenuation, angular correction, and mean chamber gain fSingleSectorFitter(0), // just for debugging // // Counters // fTotalTracks(0), // just for debugging fAcceptedTracks(0), // just for debugging fDebugCalPadRaw(0), // just for debugging fDebugCalPadCorr(0), // just for debugging fPrevIter(0) // the calibration object in its previous iteration (will not be owned by the new object, don't forget to delete it!) { // // Constructor. // G__SetCatchException(0); this->SetNameTitle(name, title); fCuts = cuts; fPrevIter = prevIter; // // Fitter initialization // fSimpleFitter = new AliTPCFitPad(8, "hyp7", ""); fSqrtFitter = new AliTPCFitPad(8, "hyp7", ""); fLogFitter = new AliTPCFitPad(8, "hyp7", ""); fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", ""); // fFitter0M = new TLinearFitter(45,"hyp44"); fFitter1M = new TLinearFitter(45,"hyp44"); fFitter2M = new TLinearFitter(45,"hyp44"); fFitter0T = new TLinearFitter(45,"hyp44"); fFitter1T = new TLinearFitter(45,"hyp44"); fFitter2T = new TLinearFitter(45,"hyp44"); // // // Add profile histograms -JUST for visualization - Not used for real calibration // // fArrayQM=new TObjArray(73); // Qmax normalized fArrayQT=new TObjArray(73); // Qtot normalized fProfileArrayQM = new TObjArray(37); // Qmax normalized versus local X fProfileArrayQT = new TObjArray(37); // Qtot normalized versus local X fProfileArrayQM2D = new TObjArray(37); // Qmax normalized versus local X and phi fProfileArrayQT2D = new TObjArray(37); // Qtot normalized versus local X and phi char hname[1000]; for (Int_t i=0; i<73; i++){ sprintf(hname,"QM_%d",i); fArrayQM->AddAt(new TH1F(hname,hname,200,0,1000),i); sprintf(hname,"QT_%d",i); fArrayQT->AddAt(new TH1F(hname,hname,200,0,1000),i); } for (Int_t i=0; i<37;i++){ sprintf(hname,"QMvsx_%d",i); fProfileArrayQM->AddAt(new TProfile(hname,hname,50,89,250),i); sprintf(hname,"QTvsx_%d",i); fProfileArrayQT->AddAt(new TProfile(hname,hname,50,89,250),i); sprintf(hname,"QM2D_%d",i); fProfileArrayQM2D->AddAt(new TProfile2D(hname,hname,50,89,250,10,-0.15,0.15),i); sprintf(hname,"QT2D_%d",i); fProfileArrayQT2D->AddAt(new TProfile2D(hname,hname,50,89,250,10,-0.15,0.15),i); } // // just for debugging -counters // fTotalTracks = 0; fAcceptedTracks = 0; fDebugCalPadRaw = new AliTPCCalPad("DebugCalPadRaw", "All clusters simply added up before correction"); fDebugCalPadCorr = new AliTPCCalPad("DebugCalPadCorr", "All clusters simply added up after correction"); // this will be gone for the a new ROOT version > v5-17-05 for (UInt_t i = 0; i < 36; i++) { fNShortClusters[i] = 0; fNMediumClusters[i] = 0; fNLongClusters[i] = 0; } } AliTPCcalibTracksGain::~AliTPCcalibTracksGain() { // // Destructor. // Info("Destructor",""); if (fSimpleFitter) delete fSimpleFitter; if (fSqrtFitter) delete fSqrtFitter; if (fLogFitter) delete fLogFitter; if (fSingleSectorFitter) delete fSingleSectorFitter; if (fDebugStream) { //fDebugStream->GetFile()->Close(); printf("Deleting debug stream object\n"); delete fDebugStream; } if (fDebugCalPadRaw) delete fDebugCalPadRaw; if (fDebugCalPadCorr) delete fDebugCalPadCorr; } void AliTPCcalibTracksGain::Terminate(){ // // Evaluate fitters and close the debug stream. // Also move or copy the debug stream, if a debugStreamPrefix is provided. // Evaluate(); if (fDebugStream) { delete fDebugStream; fDebugStream = 0; } } void AliTPCcalibTracksGain::Process(AliTPCseed* seed) { // // Main method to be called when a new seed is supposed to be processed // and be used for gain calibration. Its quality is checked before it // is added. // fTotalTracks++; if (!fCuts->AcceptTrack(seed)) return; fAcceptedTracks++; AddTrack(seed); } Long64_t AliTPCcalibTracksGain::Merge(TCollection *list) { // // Merge() merges the results of all AliTPCcalibTracksGain objects contained in // list, thus allowing a distributed computation of several files, e.g. on PROOF. // The merged results are merged with the data members of the AliTPCcalibTracksGain // object used for calling the Merge method. // The return value is 0 /*the total number of tracks used for calibration*/ if the merge // is successful, otherwise it is -1. // if (!list || list->IsEmpty()) return -1; if (!fSimpleFitter) fSimpleFitter = new AliTPCFitPad(8, "hyp7", ""); if (!fSqrtFitter) fSqrtFitter = new AliTPCFitPad(8, "hyp7", ""); if (!fLogFitter) fLogFitter = new AliTPCFitPad(8, "hyp7", ""); if (!fSingleSectorFitter) fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", ""); // just for debugging if (!fDebugCalPadRaw) fDebugCalPadRaw = new AliTPCCalPad("DebugCalPadRaw", "All clusters simply added up before correction"); if (!fDebugCalPadCorr) fDebugCalPadCorr = new AliTPCCalPad("DebugCalPadCorr", "All clusters simply added up after correction"); TIterator* iter = list->MakeIterator(); AliTPCcalibTracksGain* cal = 0; while ((cal = (AliTPCcalibTracksGain*)iter->Next())) { if (!cal->InheritsFrom(AliTPCcalibTracksGain::Class())) { Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName()); return -1; } Add(cal); } return 0; } void AliTPCcalibTracksGain::Add(AliTPCcalibTracksGain* cal) { // // Adds another AliTPCcalibTracksGain object to this object. // fSimpleFitter->Add(cal->fSimpleFitter); fSqrtFitter->Add(cal->fSqrtFitter); fLogFitter->Add(cal->fLogFitter); fSingleSectorFitter->Add(cal->fSingleSectorFitter); // // // fFitter0M->Add(cal->fFitter0M); fFitter1M->Add(cal->fFitter1M); fFitter2M->Add(cal->fFitter2M); fFitter0T->Add(cal->fFitter0T); fFitter1T->Add(cal->fFitter1T); fFitter2T->Add(cal->fFitter2T); // // // histograms // for (Int_t i=0; i<73; i++){ TH1F *his,*hism; his = (TH1F*)fArrayQM->At(i); hism = (TH1F*)cal->fArrayQM->At(i); if (his && hism) his->Add(hism); his = (TH1F*)fArrayQT->At(i); hism = (TH1F*)cal->fArrayQT->At(i); if (his && hism) his->Add(hism); } // // for (Int_t i=0; i<37; i++){ TProfile *his,*hism; his = (TProfile*)fProfileArrayQM->At(i); hism = (TProfile*)cal->fProfileArrayQM->At(i); if (his && hism) his->Add(hism); his = (TProfile*)fProfileArrayQT->At(i); hism = (TProfile*)cal->fProfileArrayQT->At(i); if (his && hism) his->Add(hism); } // // for (Int_t i=0; i<37; i++){ TProfile2D *his,*hism; his = (TProfile2D*)fProfileArrayQM2D->At(i); hism = (TProfile2D*)cal->fProfileArrayQM2D->At(i); if (his && hism) his->Add(hism); his = (TProfile2D*)fProfileArrayQT2D->At(i); hism = (TProfile2D*)cal->fProfileArrayQT2D->At(i); if (his && hism) his->Add(hism); } // // this will be gone for the a new ROOT version > v5-17-05 for (UInt_t iSegment = 0; iSegment < 36; iSegment++) { fNShortClusters[iSegment] += cal->fNShortClusters[iSegment]; fNMediumClusters[iSegment] += cal->fNMediumClusters[iSegment]; fNLongClusters[iSegment] += cal->fNLongClusters[iSegment]; } // just for debugging, remove me fTotalTracks += cal->fTotalTracks; fAcceptedTracks += cal->fAcceptedTracks; fDebugCalPadRaw->Add(cal->fDebugCalPadRaw); fDebugCalPadCorr->Add(cal->fDebugCalPadCorr); } void AliTPCcalibTracksGain::AddTrack(AliTPCseed* seed) { // // The clusters making up the track (seed) are added to various fit functions. // See AddCluster(...) for more detail. // if (!fDebugStream) fDebugStream = new TTreeSRedirector(fgkDebugStreamFileName); DumpTrack(seed); // // simple histograming part for (Int_t i=0; i<159; i++){ AliTPCclusterMI* cluster = seed->GetClusterPointer(i); if (cluster) AddCluster(cluster); } } void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster){ // // Adding cluster information to the simple histograms // No correction, fittings are applied // Float_t kThreshold=5; if (cluster->GetX()<=0) return; if (cluster->GetQ()<=kThreshold) return; // Int_t sector = cluster->GetDetector(); TH1F * his; his = GetQT(sector); if (his) his->Fill(cluster->GetQ()); his = GetQT(-1); if (his) his->Fill(cluster->GetQ()); his = GetQM(sector); if (his) his->Fill(cluster->GetMax()); his = GetQM(-1); if (his) his->Fill(cluster->GetMax()); // sector = sector%36; TProfile *prof; prof = GetProfileQT(sector); if (prof) prof->Fill(cluster->GetX(),cluster->GetQ()); prof = GetProfileQT(-1); if (prof) prof->Fill(cluster->GetX(),cluster->GetQ()); prof = GetProfileQM(sector); if (prof) prof->Fill(cluster->GetX(),cluster->GetMax()); prof = GetProfileQM(-1); if (prof) prof->Fill(cluster->GetX(),cluster->GetMax()); // Float_t phi = cluster->GetY()/cluster->GetX(); TProfile2D *prof2; prof2 = GetProfileQT2D(sector); if (prof2) prof2->Fill(cluster->GetX(),phi,cluster->GetQ()); prof2 = GetProfileQT2D(-1); if (prof2) prof2->Fill(cluster->GetX(),phi,cluster->GetQ()); prof2 = GetProfileQM2D(sector); if (prof2) prof2->Fill(cluster->GetX(),phi,cluster->GetMax()); prof2 = GetProfileQM2D(-1); if (prof2) prof2->Fill(cluster->GetX(),phi,cluster->GetMax()); // } void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster, Float_t /*momenta*/, Float_t/* mdedx*/, Int_t padType, Float_t xcenter, TVectorD& dedxQ, TVectorD& /*dedxM*/, Float_t /*fraction*/, Float_t fraction2, Float_t dedge, TVectorD& /*parY*/, TVectorD& /*parZ*/, TVectorD& meanPos) { // // Adds cluster to the appropriate fitter for later analysis. // The charge used for the fit is the maximum charge for this specific cluster or the // accumulated charge per cluster, depending on the value of fgkUseTotalCharge. // Depending on the pad size where the cluster is registered, the value will be put in // the appropriate fitter. Furthermore, for each pad size three different types of fitters // are used. The fit functions are the same for all fitters (parabolic functions), but the value // added to each fitter is different. The simple fitter gets the charge plugged in as is, the sqrt fitter // gets the square root of the charge, and the log fitter gets fgkM*(1+q/fgkM), where q is the original charge // and fgkM==25. // if (!cluster) { Error("AddCluster", "Cluster not valid."); return; } if (dedge < 3.) return; if (fraction2 > 0.7) return; //Int_t padType = GetPadType(cluster->GetX()); Double_t xx[7]; //Double_t centerPad[2] = {0}; //AliTPCFitPad::GetPadRegionCenterLocal(padType, centerPad); //xx[0] = cluster->GetX() - centerPad[0]; //xx[1] = cluster->GetY() - centerPad[1]; xx[0] = cluster->GetX() - xcenter; xx[1] = cluster->GetY(); xx[2] = xx[0] * xx[0]; xx[3] = xx[1] * xx[1]; xx[4] = xx[0] * xx[1]; xx[5] = TMath::Abs(cluster->GetZ()) - TMath::Abs(meanPos[4]); xx[6] = xx[5] * xx[5]; // // Update profile histograms // // // Update fitters // Int_t segment = cluster->GetDetector() % 36; Double_t q = fgkUseTotalCharge ? ((Double_t)(cluster->GetQ())) : ((Double_t)(cluster->GetMax())); // note: no normalization to pad size! // just for debugging Int_t row = 0; Int_t pad = 0; GetRowPad(cluster->GetX(), cluster->GetY(), row, pad); fDebugCalPadRaw->GetCalROC(cluster->GetDetector())->SetValue(row, pad, q + fDebugCalPadRaw->GetCalROC(cluster->GetDetector())->GetValue(row, pad)); // correct charge by normalising to mean charge per track q /= dedxQ[2]; // just for debugging fDebugCalPadCorr->GetCalROC(cluster->GetDetector())->SetValue(row, pad, q + fDebugCalPadCorr->GetCalROC(cluster->GetDetector())->GetValue(row, pad)); Double_t sqrtQ = TMath::Sqrt(q); Double_t logQ = fgkM * TMath::Log(1 + q / fgkM); fSimpleFitter->GetFitter(segment, padType)->AddPoint(xx, q); fSqrtFitter->GetFitter(segment, padType)->AddPoint(xx, sqrtQ); fLogFitter->GetFitter(segment, padType)->AddPoint(xx, logQ); fSingleSectorFitter->GetFitter(0, padType)->AddPoint(xx, q); // this will be gone for the a new ROOT version > v5-17-05 if (padType == kShortPads) fNShortClusters[segment]++; if (padType == kMediumPads) fNMediumClusters[segment]++; if (padType == kLongPads) fNLongClusters[segment]++; } void AliTPCcalibTracksGain::Evaluate(Bool_t robust, Double_t frac) { // // Evaluates all fitters contained in this object. // If the robust option is set to kTRUE a robust fit is performed with frac as // the minimal fraction of good points (see TLinearFitter::EvalRobust for details). // Beware: Robust fitting is much slower! // fSimpleFitter->Evaluate(robust, frac); fSqrtFitter->Evaluate(robust, frac); fLogFitter->Evaluate(robust, frac); fSingleSectorFitter->Evaluate(robust, frac); fFitter0M->Eval(); fFitter1M->Eval(); fFitter2M->Eval(); fFitter0T->Eval(); fFitter1T->Eval(); fFitter2T->Eval(); } AliTPCCalPad* AliTPCcalibTracksGain::CreateFitCalPad(UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) { // // Creates the calibration object AliTPCcalPad using fitted parameterization // TObjArray tpc(72); for (UInt_t iSector = 0; iSector < 72; iSector++) tpc.Add(CreateFitCalROC(iSector, fitType, undoTransformation, normalizeToPadSize)); return new AliTPCCalPad(&tpc); } AliTPCCalROC* AliTPCcalibTracksGain::CreateFitCalROC(UInt_t sector, UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) { // // Create the AliTPCCalROC with the values per pad // sector - sector of interest // fitType - // TVectorD par(8); if (sector < 36) { GetParameters(sector % 36, 0, fitType, par); return CreateFitCalROC(sector, 0, par, fitType, undoTransformation, normalizeToPadSize); } else { GetParameters(sector % 36, 1, fitType, par); AliTPCCalROC* roc1 = CreateFitCalROC(sector, 1, par, fitType, undoTransformation, normalizeToPadSize); GetParameters(sector % 36, 2, fitType, par); AliTPCCalROC* roc2 = CreateFitCalROC(sector, 2, par, fitType, undoTransformation, normalizeToPadSize); AliTPCCalROC* roc3 = CreateCombinedCalROC(roc1, roc2); delete roc1; delete roc2; return roc3; } } AliTPCCalROC* AliTPCcalibTracksGain::CreateFitCalROC(UInt_t sector, UInt_t padType, TVectorD &fitParam, UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) { // // This function is essentially a copy of AliTPCCalROC::CreateGlobalFitCalROC(...), with the // modifications, that the center of the region of same pad size is used as the origin // of the fit function instead of the center of the ROC. // The possibility of a linear fit is removed as well because it is not needed. // Only values for pads with the given pad size are calculated, the rest is 0. // Set undoTransformation for undoing the transformation that was applied to the // charge values before they were put into the fitter (thus allowing comparison to the original // charge values). For fitType use 0 for the simple fitter, 1 for the sqrt fitter, 2 for the log fitter. // If normalizeToPadSize is true, the values are normalized to the pad size. // Please be aware, that you even need to specify the fitType if you want to normalize to the pad size without // undoing the transformation (because normalizing involves undoing the trafo first, then normalizing, then // applying the trafo again). // Please note: The normalization to the pad size is a simple linear scaling with the pad length, which // actually doesn't describe reality! // Float_t dlx, dly; Double_t centerPad[2] = {0}; Float_t localXY[3] = {0}; AliTPCROC* tpcROC = AliTPCROC::Instance(); if ((padType == 0 && sector >= tpcROC->GetNInnerSector()) || (padType > 0 && sector < tpcROC->GetNInnerSector()) || sector >= tpcROC->GetNSector()) return 0; AliTPCCalROC* lROCfitted = new AliTPCCalROC(sector); //tpcROC->GetPositionLocal(sector, lROCfitted->GetNrows()/2, lROCfitted->GetNPads(lROCfitted->GetNrows()/2)/2, centerPad); // use this instead of the switch statement if you want to calculate the center of the ROC and not the center of the regions with the same pad size UInt_t startRow = 0; UInt_t endRow = 0; switch (padType) { case kShortPads: startRow = 0; endRow = lROCfitted->GetNrows(); break; case kMediumPads: startRow = 0; endRow = 64; break; case kLongPads: startRow = 64; endRow = lROCfitted->GetNrows(); break; } AliTPCFitPad::GetPadRegionCenterLocal(padType, centerPad); Double_t value = 0; for (UInt_t irow = startRow; irow < endRow; irow++) { for (UInt_t ipad = 0; ipad < lROCfitted->GetNPads(irow); ipad++) { tpcROC->GetPositionLocal(sector, irow, ipad, localXY); // calculate position localXY by pad and row number dlx = localXY[0] - centerPad[0]; dly = localXY[1] - centerPad[1]; value = fitParam[0] + fitParam[1]*dlx + fitParam[2]*dly + fitParam[3]*dlx*dlx + fitParam[4]*dly*dly + fitParam[5]*dlx*dly; // Let q' = value be the transformed value without any pad size corrections, // let T be the transformation and let l be the pad size // 1) don't undo transformation, don't normalize: return q' // 2) undo transformation, don't normalize: return T^{-1} q' // 3) undo transformation, normalize: return (T^{-1} q') / l // 4) don't undo transformation, normalize: return T((T^{-1} q') / l) if (!undoTransformation && !normalizeToPadSize) {/* value remains unchanged */} // (1) else { // (2), (3), (4) //calculate T^{-1} switch (fitType) { case 0: /* value remains unchanged */ break; case 1: value = value * value; break; case 2: value = (TMath::Exp(value / fgkM) - 1) * fgkM; break; default: Error("CreateFitCalROC", "Wrong fit type."); break; } if (normalizeToPadSize) value /= GetPadLength(localXY[0]); // (3) } if (!undoTransformation && normalizeToPadSize) { // (4) // calculate T switch (fitType) { case 0: /* value remains unchanged */ break; case 1: value = TMath::Sqrt(value); break; case 2: value = fgkM * TMath::Log(1 + value / fgkM); break; default: Error("CreateFitCalROC", "Wrong fit type."); break; } } lROCfitted->SetValue(irow, ipad, value); } } return lROCfitted; } AliTPCCalROC* AliTPCcalibTracksGain::CreateCombinedCalROC(const AliTPCCalROC* roc1, const AliTPCCalROC* roc2) { // // Combines the medium pad size values of roc1 with the long pad size values of roc2 into a new // AliTPCCalROC. Returns a null pointer if any one of the ROCs is an IROC; issues a warning message // if the sectors of roc1 and roc2 don't match, but still continue and use the sector of roc1 as the // sector of the new ROC. // if (!roc1 || !roc2) return 0; if ((Int_t)(roc1->GetSector()) < fgTPCparam->GetNInnerSector()) return 0; if ((Int_t)(roc2->GetSector()) < fgTPCparam->GetNInnerSector()) return 0; if (roc1->GetSector() != roc2->GetSector()) Warning("CreateCombinedCalROC", "Sector number mismatch."); AliTPCCalROC* roc = new AliTPCCalROC(roc1->GetSector()); for (UInt_t iRow = 0; iRow < 64; iRow++) { for (UInt_t iPad = 0; iPad < roc->GetNPads(iRow); iPad++) roc->SetValue(iRow, iPad, roc1->GetValue(iRow, iPad)); } for (UInt_t iRow = 64; iRow < roc->GetNrows(); iRow++) { for (UInt_t iPad = 0; iPad < roc->GetNPads(iRow); iPad++) roc->SetValue(iRow, iPad, roc2->GetValue(iRow, iPad)); } return roc; } void AliTPCcalibTracksGain::GetParameters(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitParam) { // // Puts the fit parameters for the specified segment (IROC & OROC), padType and fitType // into the fitParam TVectorD (which should contain 8 elements). // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. // Note: The fitter has to be evaluated first! // GetFitter(segment, padType, fitType)->GetParameters(fitParam); } void AliTPCcalibTracksGain::GetErrors(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitError) { // // Puts the fit parameter errors for the specified segment (IROC & OROC), padType and fitType // into the fitError TVectorD (which should contain 8 elements). // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. // Note: The fitter has to be evaluated first! // GetFitter(segment, padType, fitType)->GetErrors(fitError); fitError *= TMath::Sqrt(GetRedChi2(segment, padType, fitType)); } Double_t AliTPCcalibTracksGain::GetRedChi2(UInt_t segment, UInt_t padType, UInt_t fitType) { // // Returns the reduced chi^2 value for the specified segment, padType and fitType. // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. // Note: The fitter has to be evaluated first! // // this will be gone for the a new ROOT version > v5-17-05 Int_t lNClusters = 0; switch (padType) { case kShortPads: lNClusters = fNShortClusters[segment]; break; case kMediumPads: lNClusters = fNMediumClusters[segment]; break; case kLongPads: lNClusters = fNLongClusters[segment]; break; } return GetFitter(segment, padType, fitType)->GetChisquare()/(lNClusters - 8); } void AliTPCcalibTracksGain::GetCovarianceMatrix(UInt_t segment, UInt_t padType, UInt_t fitType, TMatrixD& covMatrix) { // // Returns the covariance matrix for the specified segment, padType, fitType. // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. // GetFitter(segment, padType, fitType)->GetCovarianceMatrix(covMatrix); } TLinearFitter* AliTPCcalibTracksGain::GetFitter(UInt_t segment, UInt_t padType, UInt_t fitType) { // // Returns the TLinearFitter object for the specified segment, padType, fitType. // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. // switch (fitType) { case kSimpleFitter: return fSimpleFitter->GetFitter(segment, padType); case kSqrtFitter: return fSqrtFitter->GetFitter(segment, padType); case kLogFitter: return fLogFitter->GetFitter(segment, padType); case 3: return fSingleSectorFitter->GetFitter(0, padType); } return 0; } Double_t AliTPCcalibTracksGain::GetPadLength(Double_t lx) { // // The function returns 0.75 for an IROC, 1. for an OROC at medium pad size position, // 1.5 for an OROC at long pad size position, -1 if out of bounds. // Double_t irocLow = fgTPCparam->GetPadRowRadiiLow(0) - fgTPCparam->GetInnerPadPitchLength()/2; Double_t irocUp = fgTPCparam->GetPadRowRadiiLow(fgTPCparam->GetNRowLow()-1) + fgTPCparam->GetInnerPadPitchLength()/2; Double_t orocLow1 = fgTPCparam->GetPadRowRadiiUp(0) - fgTPCparam->GetOuter1PadPitchLength()/2; Double_t orocUp1 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()-1) + fgTPCparam->GetOuter1PadPitchLength()/2; Double_t orocLow2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()) - fgTPCparam->GetOuter2PadPitchLength()/2; Double_t orocUp2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp()-1) + fgTPCparam->GetOuter2PadPitchLength()/2; // if IROC if (lx >= irocLow && lx <= irocUp) return 0.75; // if OROC medium pads if (lx >= orocLow1 && lx <= orocUp1) return 1.; // if OROC long pads if (lx >= orocLow2 && lx <= orocUp2) return 1.5; // if out of bounds return -1; } Int_t AliTPCcalibTracksGain::GetPadType(Double_t lx) { // // The function returns 0 for an IROC, 1 for an OROC at medium pad size position, // 2 for an OROC at long pad size position, -1 if out of bounds. // if (GetPadLength(lx) == 0.75) return 0; else if (GetPadLength(lx) == 1.) return 1; else if (GetPadLength(lx) == 1.5) return 2; return -1; } // ONLY FOR DEBUGGING PURPOSES - REMOVE ME WHEN NOT NEEDED ANYMORE Bool_t AliTPCcalibTracksGain::GetRowPad(Double_t lx, Double_t ly, Int_t& row, Int_t& pad) { // // Calculate the row and pad number when the local coordinates are given. // Returns kFALSE if the position is out of range, otherwise return kTRUE. // WARNING: This function is preliminary and probably isn't very accurate!! // Double_t irocLow = fgTPCparam->GetPadRowRadiiLow(0) - fgTPCparam->GetInnerPadPitchLength()/2; //Double_t irocUp = fgTPCparam->GetPadRowRadiiLow(fgTPCparam->GetNRowLow()-1) + fgTPCparam->GetInnerPadPitchLength()/2; Double_t orocLow1 = fgTPCparam->GetPadRowRadiiUp(0) - fgTPCparam->GetOuter1PadPitchLength()/2; //Double_t orocUp1 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()-1) + fgTPCparam->GetOuter1PadPitchLength()/2; Double_t orocLow2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()) - fgTPCparam->GetOuter2PadPitchLength()/2; //Double_t orocUp2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp()-1) + fgTPCparam->GetOuter2PadPitchLength()/2; if (GetPadType(lx) == 0) { row = (Int_t)((lx - irocLow) / fgTPCparam->GetInnerPadPitchLength()); pad = (Int_t)((ly + fgTPCparam->GetYInner(row)) / fgTPCparam->GetInnerPadPitchWidth()); } else if (GetPadType(lx) == 1) { row = (Int_t)((lx - orocLow1) / fgTPCparam->GetOuter1PadPitchLength()); pad = (Int_t)((ly + fgTPCparam->GetYOuter(row)) / fgTPCparam->GetOuterPadPitchWidth()); } else if (GetPadType(lx) == 2) { row = fgTPCparam->GetNRowUp1() + (Int_t)((lx - orocLow2) / fgTPCparam->GetOuter2PadPitchLength()); pad = (Int_t)((ly + fgTPCparam->GetYOuter(row)) / fgTPCparam->GetOuterPadPitchWidth()); } else return kFALSE; return kTRUE; } void AliTPCcalibTracksGain::DumpTrack(AliTPCseed* track) { // // Dump track information to the debug stream // Int_t rows[200]; Int_t sector[3]; Int_t npoints[3]; TVectorD dedxM[3]; TVectorD dedxQ[3]; TVectorD parY[3]; TVectorD parZ[3]; TVectorD meanPos[3]; // Int_t count=0; for (Int_t ipad = 0; ipad < 3; ipad++) { dedxM[ipad].ResizeTo(5); dedxQ[ipad].ResizeTo(5); parY[ipad].ResizeTo(3); parZ[ipad].ResizeTo(3); meanPos[ipad].ResizeTo(6); Bool_t isOK = GetDedx(track, ipad, rows, sector[ipad], npoints[ipad], dedxM[ipad], dedxQ[ipad], parY[ipad], parZ[ipad], meanPos[ipad]); if (isOK) AddTracklet(sector[ipad],ipad, dedxQ[ipad], dedxM[ipad], parY[ipad], parZ[ipad], meanPos[ipad] ); if (isOK) count++; } (*fDebugStream) << "Track" << "Track.=" << track << // track information "\n"; // // // if (count>1){ (*fDebugStream) << "TrackG" << "Track.=" << track << // track information "ncount="<GetNRowLow(); xcenter = 108.47; } if (padType == 1) { firstRow = fgTPCparam->GetNRowLow(); lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1(); xcenter = 166.60; } if (padType == 2) { firstRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1(); lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp(); xcenter = 222.6; } minRow = (lastRow - firstRow) / 2; // // Int_t nclusters = 0; Int_t nclustersNE = 0; // number of not edge clusters Int_t lastSector = -1; Float_t amplitudeQ[100]; Float_t amplitudeM[100]; Int_t rowIn[100]; Int_t index[100]; // // fitY.ClearPoints(); fitZ.ClearPoints(); for (Int_t iCluster = firstRow; iCluster < lastRow; iCluster++) { AliTPCclusterMI* cluster = track->GetClusterPointer(iCluster); if (cluster) { Int_t detector = cluster->GetDetector() ; if (lastSector == -1) lastSector = detector; if (lastSector != detector) continue; amplitudeQ[nclusters] = cluster->GetQ(); amplitudeM[nclusters] = cluster->GetMax(); rowIn[nclusters] = iCluster; nclusters++; Double_t dx = cluster->GetX() - xcenter; Double_t y = cluster->GetY(); Double_t z = cluster->GetZ(); fitY.AddPoint(&dx, y); fitZ.AddPoint(&dx, z); meanPos[0] += dx; meanPos[1] += dx; meanPos[2] += y; meanPos[3] += y*y; meanPos[4] += z; meanPos[5] += z*z; if (TMath::Abs(cluster->GetY()) < cluster->GetX()*ktany - kedgey) nclustersNE++; } } if (nclusters < minRow / 2) return kFALSE; if (nclustersNE < minRow / 2) return kFALSE; for (Int_t i = 0; i < 6; i++) meanPos[i] /= Double_t(nclusters); fitY.Eval(); fitZ.Eval(); fitY.GetParameters(parY); fitZ.GetParameters(parZ); // // calculate truncated mean // TMath::Sort(nclusters, amplitudeQ, index, kFALSE); // // // Float_t ndedx[5]; for (Int_t i = 0; i < 5; i++) { dedxQ[i] = 0; dedxM[i] = 0; ndedx[i] = 0; } // // dedx calculation // Int_t inonEdge = 0; for (Int_t i = 0; i < nclusters; i++) { Int_t rowSorted = rowIn[index[i]]; AliTPCclusterMI* cluster = track->GetClusterPointer(rowSorted); if (TMath::Abs(cluster->GetY()) > cluster->GetX()*ktany - kedgey) continue; //don't take edge clusters inonEdge++; if (inonEdge < nclustersNE * 0.5) { ndedx[0]++; dedxQ[0] += amplitudeQ[index[i]]; dedxM[0] += amplitudeM[index[i]]; } if (inonEdge < nclustersNE * 0.6) { ndedx[1]++; dedxQ[1] += amplitudeQ[index[i]]; dedxM[1] += amplitudeM[index[i]]; } if (inonEdge < nclustersNE * 0.7) { ndedx[2]++; dedxQ[2] += amplitudeQ[index[i]]; dedxM[2] += amplitudeM[index[i]]; } if (inonEdge < nclustersNE * 0.8) { ndedx[3]++; dedxQ[3] += amplitudeQ[index[i]]; dedxM[3] += amplitudeM[index[i]]; } if (inonEdge < nclustersNE * 0.9) { ndedx[4]++; dedxQ[4] += amplitudeQ[index[i]]; dedxM[4] += amplitudeM[index[i]]; } } for (Int_t i = 0; i < 5; i++) { dedxQ[i] /= ndedx[i]; dedxM[i] /= ndedx[i]; } inonEdge = 0; Float_t momenta = track->GetP(); Float_t mdedx = track->GetdEdx(); for (Int_t i = 0; i < nclusters; i++) { Int_t rowSorted = rowIn[index[i]]; AliTPCclusterMI* cluster = track->GetClusterPointer(rowSorted); if (!cluster) { printf("Problem\n"); continue; } if (TMath::Abs(cluster->GetY()) < cluster->GetX()*ktany - kedgey) inonEdge++; Float_t dedge = cluster->GetX()*ktany - TMath::Abs(cluster->GetY()); Float_t fraction = Float_t(i) / Float_t(nclusters); Float_t fraction2 = Float_t(inonEdge) / Float_t(nclustersNE); AddCluster(cluster, momenta, mdedx, padType, xcenter, dedxQ, dedxM, fraction, fraction2, dedge, parY, parZ, meanPos); (*fDebugStream) << "dEdx" << "Cl.=" << cluster << // cluster of interest "P=" << momenta << // track momenta "dedx=" << mdedx << // mean dedx - corrected for angle "IPad=" << padType << // pad type 0..2 "xc=" << xcenter << // x center of chamber "dedxQ.=" << &dedxQ << // dedxQ - total charge "dedxM.=" << &dedxM << // dedxM - maximal charge "fraction=" << fraction << // fraction - order in statistic (0,1) "fraction2=" << fraction2 << // fraction - order in statistic (0,1) "dedge=" << dedge << // distance to the edge "parY.=" << &parY << // line fit "parZ.=" << &parZ << // line fit "meanPos.=" << &meanPos << // mean position (dx, dx^2, y,y^2, z, z^2) "\n"; } (*fDebugStream) << "dEdxT" << "P=" << momenta << // track momenta "npoints="<SetAlias("dr","(250-abs(meanPos.fElements[4]))/250"); //chain->SetAlias("tz","(0+abs(parZ.fElements[1]))"); //chain->SetAlias("ty","(0+abs(parY.fElements[1]))"); //chain->SetAlias("corrg","sqrt((1+ty^2)*(1+tz^2))"); //expession fast - TString *strq0 = toolkit.FitPlane(chain,"dedxQ.fElements[2]","dr++ty++tz++dr*ty++dr*tz++ty*tz++ty^2++tz^2","IPad==0",chi2,npoints,param,covar,0,100000); Double_t xxx[100]; // // z and angular part // xxx[0] = (250.-TMath::Abs(meanPos[4]))/250.; xxx[1] = TMath::Abs(parY[1]); xxx[2] = TMath::Abs(parZ[1]); xxx[3] = xxx[0]*xxx[1]; xxx[4] = xxx[0]*xxx[2]; xxx[5] = xxx[1]*xxx[2]; xxx[6] = xxx[0]*xxx[0]; xxx[7] = xxx[1]*xxx[1]; xxx[8] = xxx[2]*xxx[2]; // // chamber part // Int_t tsector = sector%36; for (Int_t i=0;i<35;i++){ xxx[9+i]=(i==tsector)?1:0; } TLinearFitter *fitterM = fFitter0M; if (padType==1) fitterM=fFitter1M; if (padType==2) fitterM=fFitter2M; fitterM->AddPoint(xxx,dedxM[1]); // TLinearFitter *fitterT = fFitter0T; if (padType==1) fitterT = fFitter1T; if (padType==2) fitterT = fFitter2T; fitterT->AddPoint(xxx,dedxQ[1]); } TGraph *AliTPCcalibTracksGain::CreateAmpGraph(Int_t ipad, Bool_t qmax){ // // create the amplitude graph // The normalized amplitudes are extrapolated to the 0 angle (y,z) and 0 drift length // TVectorD vec; if (qmax){ if (ipad==0) fFitter0M->GetParameters(vec); if (ipad==1) fFitter1M->GetParameters(vec); if (ipad==2) fFitter2M->GetParameters(vec); }else{ if (ipad==0) fFitter0T->GetParameters(vec); if (ipad==1) fFitter1T->GetParameters(vec); if (ipad==2) fFitter2T->GetParameters(vec); } Float_t amp[36]; Float_t sec[36]; for (Int_t i=0;i<35;i++){ sec[i]=i; amp[i]=vec[10+i]+vec[0]; } amp[35]=vec[0]; Float_t mean = TMath::Mean(36,amp); for (Int_t i=0;i<36;i++){ sec[i]=i; amp[i]=(amp[i]-mean)/mean; } TGraph *gr = new TGraph(36,sec,amp); return gr; } void AliTPCcalibTracksGain::UpdateClusterParam(AliTPCClusterParam* clparam){ // // SetQ normalization parameters // // void SetQnorm(Int_t ipad, Int_t itype, TVectorD * norm); TVectorD vec; // fFitter0T->GetParameters(vec); clparam->SetQnorm(0,0,&vec); fFitter1T->GetParameters(vec); clparam->SetQnorm(1,0,&vec); fFitter2T->GetParameters(vec); clparam->SetQnorm(2,0,&vec); // fFitter0M->GetParameters(vec); clparam->SetQnorm(0,1,&vec); fFitter1M->GetParameters(vec); clparam->SetQnorm(1,1,&vec); fFitter2M->GetParameters(vec); clparam->SetQnorm(2,1,&vec); // } void AliTPCcalibTracksGain::Analyze(){ Evaluate(); }