TRD/TRDbase/AliTRDCalTrkAttach.cxx

   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 // Container for calibration parameters for AliTRDseedV1::AttachClusters()
  19 // .... Longer description of content ...
  20 //
  21 // For calibration procedure check AliTRDtrackleOflHelper::CalibrateAttach()
  22 // .... some reference to the calibration procedure .........
  23 //                                                                      //
  24 // Authors:                                                             //
  25 //   Alex Bercuci <a.bercuci@gsi.de>                                    //
  26 //                                                                      //
  27 //////////////////////////////////////////////////////////////////////////
  28
  29 #include <TFile.h>
  30 #include <TROOT.h>
  31 #include <TMath.h>
  32 #include "TH1.h"
  33 #include "TGraphErrors.h"
  34 #include "TObjArray.h"
  35
  36 #include "AliLog.h"
  37
  38 #include "AliTRDcalibDB.h"
  39 #include "AliTRDgeometry.h"
  40 #include "AliTRDCalTrkAttach.h"
  41
  42 ClassImp(AliTRDCalTrkAttach)
  43
  44 //______________________________________________________________________________
  45 AliTRDCalTrkAttach::AliTRDCalTrkAttach()
  46   :TNamed("AliTRDCalTrkAttach", "Calibration of AliTRDseedV1::AttachClusters")
  47   ,fRClikeLimit(0.65)
  48   ,fScaleCov(2.)
  49   ,fLike(NULL)
  50 {
  51 // Default constructor
  52
  53   fNsgmDy[0] = 5; fNsgmDy[1] = 7;
  54   fLikeMinRelDecrease[0] = .2; fLikeMinRelDecrease[1] = .3;
  55 }
  56
  57 //______________________________________________________________________________
  58 AliTRDCalTrkAttach::~AliTRDCalTrkAttach()
  59 {
  60 // Destructor
  61   if(fLike) delete fLike;
  62 }
  63
  64
  65 //______________________________________________________________________________
  66 Double_t AliTRDCalTrkAttach::CookLikelihood(Bool_t chg, Int_t ly, Float_t pt, Float_t phiTrk, Int_t n, Double_t dyr, Double_t dphi, Double_t sr) const
  67 {
  68 // Calculate likelihood for a segment to belong to a tracklet
  69 // Based on calibrated values
  70
  71   if(n<4){
  72     AliDebug(4, Form("Failed basic cut[s] : n[%d] ...", n));
  73     return 0.;
  74   }
  75   //check likelihood array
  76   if (!fLike || fLike->GetEntries() != AliTRDgeometry::kNlayer*kNcharge*kNcalib) {
  77     AliError("No usable AttachClusters calib object.");
  78     return 0.;
  79   }
  80   Int_t offset(ly*kNcalib*kNcharge); // offset for layer
  81   TGraphErrors *g(NULL);
  82   if(!(g = (TGraphErrors*)fLike->At(offset+Int_t(kResPos)*kNcharge+Int_t(chg)))){
  83     AliError("Failed retrieving AttachClusters graph.");
  84     return 0.;
  85   }
  86   // Interpolate p_t
  87   Int_t npts(g->GetN()), ip(0), jp(-1);
  88   Double_t x0, y0, x1, y1, dd(0.), invdx(0.), f[4]={0., 0., 0., 0.};
  89   for(Int_t kp(0); kp<npts; kp++){
  90     g->GetPoint(kp, x1, y1);
  91     if(x1>=pt){jp=kp; break;}
  92   }
  93   Bool_t boundary(kFALSE);
  94   if(jp<0){
  95     jp = npts-1;
  96     g->GetPoint(jp, x1, y1);
  97     ip = npts-1;
  98     boundary = kTRUE;
  99   }else if(jp==0){
 100     ip = 0;
 101     boundary = kTRUE;
 102   }else{
 103     ip = jp-1;
 104     g->GetPoint(ip, x0, y0);
 105     invdx = 1./(x0-x1);
 106   }
 107   // process pt dependences
 108   // +++ process dy
 109   Double_t dym = boundary?y1:((pt*(y0-y1) + (x0*y1-x1*y0))*invdx),
 110            sym = 0.5*(g->GetErrorY(ip)+g->GetErrorY(jp));
 111   dd      = (dyr - dym)/sym; dd*=dd;
 112   f[0] = TMath::Exp(-0.5*dd);
 113   // +++ process dphi
 114   if(!(g = (TGraphErrors*)fLike->At(offset+Int_t(kResAng)*kNcharge+Int_t(chg)))){
 115     AliError("Failed retrieving AttachClusters graph.");
 116     return 0.;
 117   }
 118   g->GetPoint(ip, x0, y0);g->GetPoint(jp, x1, y1);
 119   Double_t dpm = boundary?y1:((pt*(y0-y1) + (x0*y1-x1*y0))*invdx),
 120            spm = 0.5*(g->GetErrorY(ip)+g->GetErrorY(jp));
 121   dd      = (dphi - dpm)/spm; dd*=dd;
 122   f[1] = TMath::Exp(-0.5*dd);
 123   // +++ process no of clusters
 124   if(!(g = (TGraphErrors*)fLike->At(offset+Int_t(kNclMean)*kNcharge+Int_t(chg)))){
 125     AliError("Failed retrieving AttachClusters graph.");
 126     return 0.;
 127   }
 128   g->GetPoint(ip, x0, y0);g->GetPoint(jp, x1, y1);
 129   Double_t nm = boundary?y1:((pt*(y0-y1) + (x0*y1-x1*y0))*invdx);
 130   f[2] = (nm-TMath::Abs(n-nm))/nm;
 131
 132   // process phi dependences
 133   // +++ process <s>/s
 134   if(!(g = (TGraphErrors*)fLike->At(offset+Int_t(kSigma)*kNcharge+Int_t(chg)))){
 135     AliError("Failed retrieving AttachClusters graph.");
 136     return 0.;
 137   }
 138   // Interpolate phi [deg]
 139   npts=g->GetN(); jp=-1;
 140   for(Int_t kp(0); kp<npts; kp++){
 141     g->GetPoint(kp, x1, y1);
 142     if(x1>=phiTrk){jp=kp; break;}
 143   }
 144   if(jp<0){
 145     jp = npts-1;
 146     g->GetPoint(jp, x1, y1);
 147     ip = jp;
 148     boundary = kTRUE;
 149   }else if(jp==0){
 150     ip = jp;
 151     boundary = kTRUE;
 152   }else{
 153     ip = jp-1;
 154     g->GetPoint(ip, x0, y0);
 155     invdx = 1./(x0-x1);
 156     boundary = kFALSE;
 157   }
 158   Double_t sm = boundary?y1:((phiTrk*(y0-y1) + (x0*y1-x1*y0))*invdx),
 159            ssm = 0.5*(g->GetErrorY(ip)+g->GetErrorY(jp));
 160   dd      = (sr - sm)/ssm; dd*=dd;
 161   f[3] = 1.;//TMath::Exp(-0.5*dd);
 162
 163   // Calculate likelihood
 164   Double_t length = f[0]*f[0]+f[1]*f[1]+f[2]*f[2]+f[3]*f[3];
 165   length = TMath::Sqrt(length);
 166   Double_t cosTht = f[0]+f[1]+f[2]+f[3];
 167   cosTht /= (4.*length);
 168   AliDebug(2, Form("Like[%5.3f] ThtLike[%6.2f](deg)\n"
 169     "    F_dy (%+6.2f)=%4.2f\n"
 170     "    F_phi(%+6.2f)=%4.2f\n"
 171     "    F_ncl(%+6d)=%4.2f\n"
 172     "    F_<s>(%+6.2f)=%4.2f",
 173     length, TMath::ACos(cosTht)*TMath::RadToDeg(),
 174     dyr, f[0], dphi, f[1], n, f[2], sr, f[3]));
 175
 176   return length;
 177 }
 178
 179 //______________________________________________________________________________
 180 void AliTRDCalTrkAttach::Help()
 181 {
 182 // Display help message
 183   printf(
 184     "Draw likelihood distribution. Possible options are of the form \"lt\"\n"
 185     "where \"l\" is the layer number [0-5] and\n"
 186     "      \"t\" is the type. This is one of the following\n"
 187     "            \"y\" - r-phi roads\n"
 188     "            \"p\" - angular (deflection) roads\n"
 189     "            \"n\" - number of clusters\n"
 190     "            \"s\" - cluster spread\n");
 191 }
 192
 193 //______________________________________________________________________________
 194 void AliTRDCalTrkAttach::Draw(Option_t* opt)
 195 {
 196 // Draw likelihood distribution. Possible options are of the form "lt"
 197 // where "l" is the layer number [0-5] and
 198 //       "t" is the type. This is one of the following
 199 //             "y" - r-phi roads
 200 //             "p" - angular (deflection) roads
 201 //             "n" - number of clusters
 202 //             "s" - cluster spread
 203   if (!fLike || fLike->GetEntries() != AliTRDgeometry::kNlayer*kNcharge*kNcalib) {
 204     AliError("No likelihood distributions stored");
 205     return;
 206   }
 207   if(!opt){
 208     Help();
 209     return;
 210   }
 211   Int_t ly(-1); Char_t cly[] = {'0','1','2','3','4','5'};
 212   for(Int_t ily(0); ily<AliTRDgeometry::kNlayer; ily++){
 213     if(opt[0] == cly[ily]){
 214       ly = ily; break;
 215     }
 216   }
 217   if(ly<0){
 218     Help();
 219     return;
 220   }
 221   Int_t typ(-1); const Char_t ctyp[] = {'y', 'p', 's', 'n'};
 222   for(Int_t it(0); it<kNcalib; it++){
 223     if(opt[1] == ctyp[it]){
 224       typ = it; break;
 225     }
 226   }
 227   if(typ<0){
 228     Help();
 229     return;
 230   }
 231   Int_t offset(ly*kNcalib*kNcharge+typ*kNcharge); // offset for layer and type
 232   TGraphErrors *g[2]={(TGraphErrors*)fLike->At(offset), (TGraphErrors*)fLike->At(offset+1)};
 233   if(!g[0] || !g[1]){
 234     AliError(Form("Failed retrieving graphs for Ly[%d] Typ[%c]", ly, ctyp[typ]));
 235     return;
 236   }
 237
 238   // draw !!
 239   const Char_t *ttyp[] = {"#Deltay [cm]", "#Delta#phi [deg]", "#Delta#sigma/<#sigma>", "n_{cl}^{tracklet}"};
 240   g[0]->Draw("apl"); g[0]->SetLineColor(kBlue); g[0]->SetMarkerColor(kBlue);
 241   g[1]->Draw("pl");  g[1]->SetLineColor(kRed); g[1]->SetMarkerColor(kRed);
 242   TH1 *h=g[0]->GetHistogram();
 243   h->GetYaxis()->SetTitle(ttyp[typ]);
 244   h->GetYaxis()->SetLabelFont(52);
 245   h->GetYaxis()->SetTitleFont(62);
 246   h->GetXaxis()->SetTitle("p_{t} [GeV/c]");
 247   h->GetXaxis()->SetLabelFont(52);
 248   h->GetXaxis()->SetTitleFont(62);
 249 }
 250
 251 //______________________________________________________________________________
 252 Bool_t AliTRDCalTrkAttach::LoadReferences(const Char_t *file)
 253 {
 254 // Load calibration data from file
 255
 256   if(!file || !TFile::Open(file)){
 257     AliError("Parametrization file missing or unreadable.");
 258     return kFALSE;
 259   }
 260   TGraphErrors *g(NULL);
 261   Char_t co[kNcalib] = {'y', 'p', 's', 'n'},
 262          cs[2] = {'n', 'p'};
 263
 264   if(fLike) fLike->Clear();
 265   else fLike = new TObjArray(AliTRDgeometry::kNlayer*kNcharge*kNcalib);
 266   for(Int_t ily(0); ily<AliTRDgeometry::kNlayer; ily++){
 267     for(Int_t icalib(0); icalib<kNcalib; icalib++){
 268       for(Int_t isgn(0); isgn<kNcharge; isgn++){
 269         if(!(g = (TGraphErrors*)gFile->Get(Form("%c%c%d", co[icalib], cs[isgn], ily)))) return kFALSE;
 270         fLike->AddAt(g, kNcharge*(ily*kNcalib+icalib)+isgn);
 271       }
 272     }
 273   }
 274   return kTRUE;
 275 }