T0/AliT0CalibWalk.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 /* $Id$ */
  17
  18 ///////////////////////////////////////////////////////////////////////////////
  19 //                                                                           //
  20 // class T0 walk correction      Alla Maevskaya alla@inr.ru 20.11.2007        //
  21 //                                                                           //
  22 ///////////////////////////////////////////////////////////////////////////////
  23
  24 #include "AliT0CalibWalk.h"
  25 #include "AliLog.h"
  26
  27 #include <TObjArray.h>
  28 #include <TGraph.h>
  29 #include <TFile.h>
  30 #include <TH2F.h>
  31 #include <TMath.h>
  32 #include <TSystem.h>
  33 #include <Riostream.h>
  34 #include <TSpectrum.h>
  35 #include <TProfile.h>
  36
  37
  38 ClassImp(AliT0CalibWalk)
  39
  40 //________________________________________________________________
  41   AliT0CalibWalk::AliT0CalibWalk():   TNamed(),
  42                                       fWalk(0),
  43                                       fAmpLEDRec(0)
  44 {
  45   //
  46 }
  47
  48 //________________________________________________________________
  49 AliT0CalibWalk::AliT0CalibWalk(const char* name):TNamed(),
  50                                       fWalk(0),
  51                                       fAmpLEDRec(0)
  52 {
  53   TString namst = "Calib_";
  54   namst += name;
  55   SetName(namst.Data());
  56   SetTitle(namst.Data());
  57
  58 }
  59
  60 //________________________________________________________________
  61 AliT0CalibWalk::AliT0CalibWalk(const AliT0CalibWalk& calibda) :
  62   TNamed(calibda),
  63   fWalk(0),
  64   fAmpLEDRec(0)
  65
  66 {
  67 // copy constructor
  68   SetName(calibda.GetName());
  69   SetTitle(calibda.GetName());
  70
  71
  72 }
  73
  74 //________________________________________________________________
  75 AliT0CalibWalk &AliT0CalibWalk::operator =(const AliT0CalibWalk& calibda)
  76 {
  77 // assignment operator
  78   SetName(calibda.GetName());
  79   SetTitle(calibda.GetName());
  80
  81   return *this;
  82 }
  83
  84 //________________________________________________________________
  85 AliT0CalibWalk::~AliT0CalibWalk()
  86 {
  87   //
  88 }
  89 //________________________________________________________________
  90 void AliT0CalibWalk::SetWalk(Int_t ipmt)
  91 {
  92   //read QTC walk graph from external file
  93
  94   Int_t mv, ps;
  95   Int_t x[70000], y[70000], index[70000];
  96   Float_t time[10000],amplitude[10000];
  97   string buffer;
  98   Bool_t down=false;
  99
 100   const char * filename = gSystem->ExpandPathName("$ALICE_ROOT/T0/data/CFD-Amp.txt");
 101   ifstream inFile(filename);
 102   if(!inFile) AliError(Form("Cannot open file %s !",filename));
 103
 104   Int_t i=0;
 105   while(getline(inFile,buffer)){
 106     inFile >> ps >> mv;
 107
 108     x[i]=ps; y[i]=mv;
 109     i++;
 110   }
 111   inFile.close();
 112   cout<<" number of data "<<i<<endl;
 113
 114   TMath::Sort(i, y, index,down);
 115   Int_t amp=0, iin=0, isum=0, sum=0;
 116   Int_t ind=0;
 117   for (Int_t ii=0; ii<i; ii++)
 118     {
 119       ind=index[ii];
 120       if(y[ind] == amp)
 121         {
 122           sum +=x[ind];
 123           iin++;
 124         }
 125       else
 126         {
 127           if(iin>0)
 128             time[isum] = Float_t (sum/(iin));
 129           else
 130             time[isum] =Float_t (x[ind]);
 131           amplitude[isum] = Float_t (amp);
 132           amp=y[ind];
 133           iin=0;
 134           isum++;
 135           sum=0;
 136         }
 137     }
 138
 139   inFile.close();
 140
 141   TGraph* gr = new TGraph(isum, amplitude, time);
 142   fWalk.AddAtAndExpand(gr,ipmt);
 143 }
 144
 145 //________________________________________________________________
 146
 147 void AliT0CalibWalk::SetAmpLEDRec(Int_t ipmt)
 148 {
 149   // read LED walk from external file
 150   Float_t mv, ps;
 151   Float_t x[100], y[100];
 152   string buffer;
 153
 154   const char * filename = gSystem->ExpandPathName("$ALICE_ROOT/T0/data/CFD-LED.txt");
 155    ifstream inFile(filename);
 156   if(!inFile) {AliError(Form("Cannot open file %s !",filename));}
 157
 158   inFile >> mv>>ps;
 159   Int_t i=0;
 160
 161   while(getline(inFile,buffer)){
 162     x[i]=mv; y[i]=ps;
 163     inFile >> mv >> ps;
 164     i++;
 165   }
 166   inFile.close();
 167   Float_t y1[100], x1[100];
 168   for (Int_t ir=0; ir<i; ir++){
 169     y1[ir]=y[i-ir]; x1[ir]=x[i-ir];}
 170   TGraph* gr = new TGraph(i,y1,x1);
 171   fAmpLEDRec.AddAtAndExpand(gr,ipmt);
 172
 173 }
 174
 175 //________________________________________________________________
 176
 177 void AliT0CalibWalk::MakeWalkCorrGraph(const char *laserFile)
 178 {
 179   //make walk corerction for preprocessor
 180
 181   TFile *gFile = TFile::Open(laserFile);
 182
 183   Int_t npeaks = 20;
 184   Double_t sigma=3.;
 185   Bool_t down = false;
 186
 187   Int_t index[20];
 188   Char_t buf1[20], buf2[20],buf3[20];
 189
 190
 191   for (Int_t i=0; i<24; i++)
 192     {
 193       if(i>11)  {  sprintf(buf1,"T0_A_%i_CFD",i+1-12); }
 194
 195       if(i<12)  { sprintf(buf1,"T0_C_%i_CFD",i+1); }
 196
 197       sprintf(buf2,"CFDvsQTC%i",i+1);
 198       sprintf(buf3,"CFDvsLED%i",i+1);
 199       // cout<<buf1<<" "<<buf2<<" "<<buf3<<endl;
 200       TH2F *qtcVScfd = (TH2F*) gFile->Get(buf2);
 201       //        qtcVScfd->Print();
 202       TH2F *ledVScfd = (TH2F*) gFile->Get(buf3);
 203       TH1F *cfd = (TH1F*) gFile->Get(buf1);
 204
 205       TSpectrum *s = new TSpectrum(2*npeaks,1);
 206       Int_t nfound = s->Search(cfd,sigma,"goff",0.1);
 207       //      cout<<"Found "<<nfound<<" peaks sigma "<<sigma<<endl;;
 208        if(nfound!=0){
 209         Float_t *xpeak = s->GetPositionX();
 210         TMath::Sort(nfound, xpeak, index,down);
 211         Float_t xp = xpeak[index[0]];
 212         Int_t xbin = cfd->GetXaxis()->FindBin(xp);
 213         Float_t yp = cfd->GetBinContent(xbin);
 214         //      cout<<"xbin = "<<xbin<<"\txpeak = "<<xpeak[1]<<"\typeak = "<<yp<<endl;
 215         Float_t hmax = xp+20*sigma;
 216         Float_t hmin = xp-20*sigma;
 217       }
 218         //QTC
 219         Int_t nbins= qtcVScfd->GetXaxis()->GetNbins();
 220         TProfile *prY = qtcVScfd->ProfileX();
 221         Float_t hmin = qtcVScfd->GetYaxis()->GetXmin();
 222         Float_t hmax = qtcVScfd->GetYaxis()->GetXmax();
 223         prY->SetMaximum(hmax);
 224         prY->SetMinimum(hmin);
 225         Int_t np=0; //=nbins/5;
 226         Double_t *xx = new Double_t[nbins];
 227         Double_t *yy = new Double_t[nbins];
 228         Int_t ng=0;
 229         Double_t yg=0;
 230         for (Int_t ip=1; ip<nbins; ip++)
 231           {
 232
 233             if(ip%5 != 0 )
 234               {
 235               //              if (prY->GetBinContent(ip) !=0)
 236                 if (prY->GetBinContent(ip) >hmin && prY->GetBinContent(ip)<hmax){
 237                 yg +=prY->GetBinContent(ip);
 238                 ng++;
 239                 //      cout<<" ip "<<ip<<" Y "<<yg<<endl;
 240                           }
 241             }
 242             else {
 243               if (/*prY->GetBinContent(ip) >hmin && prY->GetBinContent(ip)<hmax &&*/ ng>0){
 244                 //              xx[ip/5] = Float_t (prY->GetBinCenter(ip));
 245                 //      yy[ip/5] = yg/ng;
 246                 xx[np] = Float_t (prY->GetBinCenter(ip));
 247                 yy[np] = yg/ng;
 248                 yg=0;
 249                 ng=0;
 250                 //      cout<<" ip "<<ip<<" np "<<np<<" XX "<<xx[np]<<" YY "<<yy[np]<<endl;
 251                 np++;
 252               }
 253             }
 254           }
 255
 256         TGraph *gr = new TGraph(np,xx,yy);
 257         gr->SetMinimum(hmin);
 258         gr->SetMaximum(hmax);
 259         fWalk.AddAtAndExpand(gr,i);
 260
 261         //LED
 262         Int_t nbinsled= ledVScfd->GetXaxis()->GetNbins();
 263         cout<<" nbins led "<<nbinsled<<endl;
 264         TProfile *prledY = ledVScfd->ProfileX();
 265
 266         Int_t npled = 0; //nbinsled/5;
 267         Double_t *xxled = new Double_t[nbinsled];
 268         Double_t *yyled = new Double_t[nbinsled];
 269         Int_t ngled=0;
 270         Double_t ygled=0;
 271         for (Int_t ip=1; ip<nbinsled; ip++)
 272           {
 273             if(ip%5 != 0 ) {
 274               if (prledY->GetBinContent(ip) !=0)
 275                 ygled +=prledY->GetBinContent(ip);
 276               ngled++;}
 277             else {
 278               xxled[npled] = Float_t (prledY->GetBinCenter(ip));
 279               yyled[npled] = ygled/ngled;
 280               ygled=0;
 281               ngled=0;
 282               npled++;
 283             }
 284           }
 285
 286         TGraph *grled = new TGraph(npled,xxled,yyled);
 287         grled->SetMinimum(hmin);
 288         grled->SetMaximum(hmax);
 289         fAmpLEDRec.AddAtAndExpand(grled,i);
 290
 291         delete [] xx;
 292         delete [] yy;
 293         delete [] xxled;
 294         delete [] yyled;
 295         delete prY;
 296         delete prledY;
 297
 298         // }
 299       delete cfd;
 300       delete qtcVScfd;
 301       delete ledVScfd;
 302     }
 303
 304 }
 305