]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/AliMUONGain.cxx
Adding the ShowConfig function
[u/mrichter/AliRoot.git] / MUON / AliMUONGain.cxx
index 5244914c9e474d0a9c3de4f8706ec06390c126fc..b42e01dbed2af038f8c2b92754618a5dd79a445f 100644 (file)
 
 // functions
 
-
-//______________________________________________________________________________
-Double_t funcLin (const Double_t *x, const Double_t *par)
-{
-  /// Linear function
-  return par[0] + par[1]*x[0];
-}
-
-//______________________________________________________________________________
-Double_t funcParabolic (const Double_t *x, const Double_t *par)
-{
-  /// Parabolic function
-  return par[0]*x[0]*x[0];
-}
-
-//______________________________________________________________________________
-Double_t funcCalib (const Double_t *x, const Double_t *par)  
-{
-  /// Calibration function
-  Double_t xLim= par[3];
-
-  if(x[0] <= xLim) return par[0] + par[1]*x[0];
-
-  Double_t yLim = par[0]+ par[1]*xLim;
-  return yLim + par[1]*(x[0] - xLim) + par[2]*(x[0] - xLim)*(x[0] - xLim);
+namespace {
+  
+  //______________________________________________________________________________
+  Double_t funcLin (const Double_t *x, const Double_t *par)
+  {
+    /// Linear function
+    return par[0] + par[1]*x[0];
+  }
+  
+  //______________________________________________________________________________
+  Double_t funcParabolic (const Double_t *x, const Double_t *par)
+  {
+    /// Parabolic function
+    return par[0]*x[0]*x[0];
+  }
+  
+  //______________________________________________________________________________
+  Double_t funcCalib (const Double_t *x, const Double_t *par)  
+  {
+    /// Calibration function
+    Double_t xLim= par[3];
+    
+    if(x[0] <= xLim) return par[0] + par[1]*x[0];
+    
+    Double_t yLim = par[0]+ par[1]*xLim;
+    return yLim + par[1]*(x[0] - xLim) + par[2]*(x[0] - xLim)*(x[0] - xLim);
+  }
+  
 }
 
-
 /// \cond CLASSIMP
 ClassImp(AliMUONGain)
 /// \endcond
@@ -89,6 +91,7 @@ ClassImp(AliMUONGain)
 AliMUONGain::AliMUONGain()
 : AliMUONPedestal(),
 fInjCharge(0), 
+fRootDataFileName(),
 fnInit(1),
 fnEntries(11),
 fnbpf1(6),
@@ -97,9 +100,7 @@ fPlotLevel(0)
 {
 /// Default constructor
 
-  sprintf(fRootDataFileName," "); //Gain
 }
-//  AliMUONPedestal& operator=(const AliMUONPedestal& other); Copy ctor
   
 //______________________________________________________________________________
 AliMUONGain::~AliMUONGain()
@@ -115,7 +116,7 @@ TString AliMUONGain::WriteDummyHeader(void)
   ostringstream stream;
   stream<<"//DUMMY FILE (to prevent Shuttle failure)"<< endl;
   stream<<"//================================================" << endl;
-  stream<<"//       MUONTRKda: Calibration run  " << endl;
+  stream<<"//       MUONTRKGAINda: Calibration run  " << endl;
   stream<<"//================================================" << endl;
   stream<<"//   * Run           : " << fRunNumber << endl; 
   stream<<"//   * Date          : " << fDate->AsString("l") <<endl;
@@ -128,11 +129,11 @@ TString AliMUONGain::WriteDummyHeader(void)
 //______________________________________________________________________________
 void AliMUONGain::MakePedStoreForGain(TString shuttleFile)
 {
-///
+/// Store Pedmean and sigma to pedestal-like ascii file
 
   ofstream fileout;
   TString tempstring;
-  Char_t flatFile[256]="";
+  TString flatFile;
   TString outputFile;
   
   // Store pedestal map in root file
@@ -141,15 +142,12 @@ void AliMUONGain::MakePedStoreForGain(TString shuttleFile)
   // write dummy ascii file -> Shuttle
   if(fIndex<fnEntries)
     {  
-      //     fileout.open(shuttleFile.Data());
-      //     tempstring = WriteDummyHeader();
-      //     fileout << tempstring;
       FILE *pfilew=0;
       pfilew = fopen (shuttleFile,"w");
 
       fprintf(pfilew,"//DUMMY FILE (to prevent Shuttle failure)\n");
       fprintf(pfilew,"//================================================\n");
-      fprintf(pfilew,"//       MUONTRKda: Calibration run  \n");
+      fprintf(pfilew,"//       MUONTRKGAINda: Calibration run  \n");
       fprintf(pfilew,"//=================================================\n");
       fprintf(pfilew,"//   * Run           : %d \n",fRunNumber); 
       fprintf(pfilew,"//   * Date          : %s \n",fDate->AsString("l"));
@@ -165,9 +163,9 @@ void AliMUONGain::MakePedStoreForGain(TString shuttleFile)
   if(fPrintLevel>0)
     {
       // compute and store mean DAC values (like pedestals)
-      sprintf(flatFile,"%s.ped",fprefixDA);
+      flatFile = Form("%s.ped",fPrefixDA.Data());
       outputFile=flatFile;
-      cout << "\n" << fprefixDA << " : Flat file  generated  : " << flatFile << "\n";
+      cout << "\n" << fPrefixDA.Data() << " : Flat file  generated  : " << flatFile.Data() << "\n";
       if (!outputFile.IsNull())  
       {
         ofstream out(outputFile.Data());
@@ -181,7 +179,7 @@ void AliMUONGain::MakePedStoreForGain(TString shuttleFile)
   if (fIndex==1) {
     mode = "RECREATE";
   }
-  TFile* histoFile = new TFile(fRootDataFileName, mode.Data(), "MUON Tracking Gains");
+  TFile* histoFile = new TFile(fRootDataFileName.Data(), mode.Data(), "MUON Tracking Gains");
 
   // second argument should be the injected charge, taken from config crocus file
   // put also info about run number could be usefull
@@ -211,25 +209,29 @@ void AliMUONGain::MakePedStoreForGain(TString shuttleFile)
 //______________________________________________________________________________
 TString AliMUONGain::WriteGainHeader(Int_t nInit, Int_t nEntries, Int_t nbpf2, Int_t *numrun, Double_t *injCharge) 
 {
-///
+/// Header of the calibration output file
 
   ostringstream stream;
 
 
-  stream<<"//================================================" << endl;
-  stream<<"//  Calibration file calculated by MUONTRKda" << endl;
-  stream<<"//=================================================" << endl;
+  stream<<"//=======================================================" << endl;
+  stream<<"//      Calibration file calculated by " << fPrefixDA.Data() <<endl;
+  stream<<"//=======================================================" << endl;
   stream<<"//   * Run           : " << fRunNumber << endl; 
   stream<<"//   * Date          : " << fDate->AsString("l") <<endl;
   stream<<"//   * Statictics    : " << fNEvents << endl;
-  stream<<"//   * # of MANUS    : " << fNManu << endl;
+  if(fConfig)
+  stream<<"//   * # of MANUS    : " << fNManuConfig << " read in the Det. config. " << endl;
+  stream<<"//   * # of MANUS    : " << fNManu << " read in raw data " << endl;
+  stream<<"//   * # of MANUS    : " << fNChannel/64 << " written in calibration file " << endl;
   stream<<"//   * # of channels : " << fNChannel << endl;
-  stream<<"//-------------------------------------------------" << endl;
+  stream<<"//-------------------------------------------------------" << endl;
 
   if(nInit==0)
     stream<<"//  "<< nEntries <<" DAC values  fit: "<< fnbpf1 << " pts (1st order) " << nbpf2 << " pts (2nd order)" << endl;
   if(nInit==1)
     stream<<"//  "<< nEntries <<" DAC values  fit: "<< fnbpf1 << " pts (1st order) " << nbpf2 << " pts (2nd order) DAC=0 excluded" << endl;
+  stream<<"//   *  nInit = " << nInit << "  *  f1nbp = " << fnbpf1 << "  *  f2nbp = " <<  nbpf2 << endl; 
 
   stream<<"//   RUN     DAC   " << endl;
   stream<<"//-----------------" << endl;
@@ -246,7 +248,7 @@ TString AliMUONGain::WriteGainHeader(Int_t nInit, Int_t nEntries, Int_t nbpf2, I
 //______________________________________________________________________________
 TString AliMUONGain::WriteGainData(Int_t BP, Int_t Manu, Int_t ch, Double_t p1, Double_t p2, Int_t threshold, Int_t q) 
 {
-///
+/// Write calibration parameters per channel
 
   ostringstream stream("");
   stream << Form("%4i %5i %2i %7.4f %10.3e %4i %2x\n",BP,Manu,ch,p1,p2,threshold,q);
@@ -261,14 +263,14 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   ofstream fileout;
   ofstream filcouc;
   TString tempstring;  
-  Char_t filename[256]
+  TString filename
 
   Double_t goodA1Min =  0.5;
   Double_t goodA1Max =  2.;
-  //     Double_t goodA1Min =  0.7;
-  //     Double_t goodA1Max =  1.7;
-  Double_t goodA2Min = -0.5E-03;
-  Double_t goodA2Max =  1.E-03;
+//   Double_t goodA2Min = -0.5E-03;
+//   Double_t goodA2Max =  1.E-03;
+  Double_t goodA2Min = -0.5E-01; // changed 28/10/2009 (JLC) <=> enlarged condition on a2
+  Double_t goodA2Max =  1.E-01;
   // Table for uncalibrated  buspatches and manus
   THashList* uncalBuspatchManuTable = new THashList(1000,2);
 
@@ -281,7 +283,8 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   // Fit with a polynomial fct
   // store the result in a flat file.
 
-  TFile*  histoFile = new TFile(fRootDataFileName);
+  if(fIndex==0)cout << " Root data file = " << fRootDataFileName.Data() << endl;  
+  TFile*  histoFile = new TFile(fRootDataFileName.Data());
 
   AliMUON2DMap* map[11];
   AliMUONVCalibParam* ped[11];
@@ -314,7 +317,7 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   for ( Int_t i=0 ; i<11 ; i++) {injCharge[i]=0.;injChargeErr[i]=1.;};
 
   // some print
-  cout<<"\n ********  MUONTRKda for Gain computing (Last Run = " << fRunNumber << ") ********\n" << endl;
+  cout<<"\n ********  MUONTRKGAINda for Gain computing (Last Run = " << fRunNumber << ") ********\n" << endl;
   cout<<" * Date          : " << fDate->AsString("l") << "\n" << endl;
   cout << " Entries = " << nEntries << " DAC values \n" << endl; 
   for (Int_t i = 0; i < nEntries; ++i) {
@@ -329,8 +332,8 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   //  print out in .log file
 
   (*fFilcout)<<"\n\n//=================================================" << endl;
-  (*fFilcout)<<"//    MUONTRKda: Gain Computing  Run = " << fRunNumber << endl;
-  (*fFilcout)<<"//    RootDataFile  = "<< fRootDataFileName << endl;
+  (*fFilcout)<<"//    MUONTRKGAINda: Gain Computing  Run = " << fRunNumber << endl;
+  (*fFilcout)<<"//    RootDataFile  = "<< fRootDataFileName.Data() << endl;
   (*fFilcout)<<"//=================================================" << endl;
   (*fFilcout)<<"//* Date          : " << fDate->AsString("l") << "\n" << endl;
 
@@ -340,9 +343,9 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
     FILE *pfilen = 0;
     if(fPrintLevel>1)
       {
-        sprintf(filename,"%s.param",fprefixDA);
-        cout << " Second fit parameter file        = " << filename << "\n";
-        pfilen = fopen (filename,"w");
+        filename=Form("%s.param",fPrefixDA.Data());
+        cout << " Second fit parameter file        = " << filename.Data() << "\n";
+        pfilen = fopen (filename.Data(),"w");
 
         fprintf(pfilen,"//===================================================================\n");
         fprintf(pfilen,"//  BP MANU CH. par[0]     [1]     [2]     [3]      xlim          P(chi2) p1        P(chi2)2  p2\n");
@@ -363,7 +366,7 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   FILE *pfilep = 0;
   if(fPrintLevel>1)
     {
-      sprintf(filename,"%s.peak",fprefixDA);
+      filename=Form("%s.peak",fPrefixDA.Data());
       cout << " File containing Peak mean values = " << filename << "\n";
       pfilep = fopen (filename,"w");
 
@@ -389,7 +392,7 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   Int_t q = 0;
   Int_t p1 =0;
   Int_t p2 =0;
-  Double_t gain=0; 
+  Double_t gain=10.; // max value (= bad value)
   Double_t capa=0.2; // internal capacitor (pF)
 
   //  plot out 
@@ -398,8 +401,8 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   TTree* tg = 0x0;
   if(fPlotLevel>0)
     {
-      sprintf(fHistoFileName,"%s.root",fprefixDA);
-     gainFile = new TFile(fHistoFileName,"RECREATE","MUON Tracking gains");
+      fHistoFileName=Form("%s.root",fPrefixDA.Data());
+      gainFile = new TFile(fHistoFileName.Data(),"RECREATE","MUON Tracking gains");
       tg = new TTree("tg","TTree avec class Manu_DiMu");
 
       tg->Branch("bp",&busPatchId, "busPatchId/I");
@@ -438,10 +441,12 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   Double_t xp[11], xpErr[11], yp[11], ypErr[11];
 
   Int_t uncalcountertotal=0 ;
+  Int_t unparabolicfit=0;
 
   while ( ( p = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
     {
       ped[0]  = p;
+      unparabolicfit=0;
 
       busPatchId = p->ID0();
       manuId     = p->ID1();
@@ -468,7 +473,6 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
            n++;
          }
 
-
          // print_peak_mean_values
          if(fPrintLevel>1)
            {
@@ -486,16 +490,24 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
 
          // 1. - linear fit over gAlinbpf1 points
 
-         Double_t par[4] = {0.,0.5,0.,kADCMax};
+         Double_t par[4] = {0.,0.5,0.,ADCMax()};
          Int_t nbs   = nEntries - fnInit;
          if(nbs < fnbpf1)fnbpf1=nbs;
 
          Int_t fitproceed=1;
+         Int_t nbpf2Dynamic=nbpf2;
+         Int_t ADClimit=4090; // when RMS < 0.5 (in other cases mean values forced to 4095, see DA_PED)
          for (Int_t j = 0; j < nbs; ++j)
            {
              Int_t k = j + fnInit;
              x[j]    = pedMean[k];
-             if(x[j]==0. || x[j]== kADCMax)fitproceed=0;
+             if(x[j]<=0.){fitproceed=0; break;}
+             //              if(x[j]>= ADCMax())
+             if(x[j]>= ADClimit)
+               {
+                 if(j < nbs-1){fitproceed=0; break;}
+                 else { nbpf2Dynamic=nbpf2-1; break;}
+               }
              xErr[j] = pedSigma[k];
              y[j]    = injCharge[k];
              yErr[j] = injChargeErr[k];
@@ -508,7 +520,7 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
          if(fitproceed)
            {
                      
-             TF1 *f1 = new TF1("f1",funcLin,0.,kADCMax,2);
+             TF1 *f1 = new TF1("f1",funcLin,0.,ADCMax(),2);
              graphErr = new TGraphErrors(fnbpf1, x, y, xErr, yErr);
 
              f1->SetParameters(0,0);
@@ -531,10 +543,10 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
              a1 = par[1];
 
              // 2. - Translation : new origin (xLim, yLim) + parabolic fit over nbf2 points
-
-             if(nbpf2 > 1)
+             //checking:         if(busPatchId ==1841 && manuId==4)nbpf2Dynamic=2;
+             if(nbpf2Dynamic > 2)
                {
-                 for (Int_t j = 0; j < nbpf2; j++)
+                 for (Int_t j = 0; j < nbpf2Dynamic; j++)
                    {
                      Int_t k  = j + (fnInit + fnbpf1) - 1;
                      xp[j]    = pedMean[k] - xLim;
@@ -544,8 +556,8 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
                      ypErr[j] = injChargeErr[k];
                    }
 
-                 TF1 *f2 = new TF1("f2",funcParabolic,0.,kADCMax,1);
-                 graphErr = new TGraphErrors(nbpf2, xp, yp, xpErr, ypErr);
+                 TF1 *f2 = new TF1("f2",funcParabolic,0.,ADCMax(),1);
+                 graphErr = new TGraphErrors(nbpf2Dynamic, xp, yp, xpErr, ypErr);
 
                  graphErr->Fit(f2,"RQ");
                  chi2P2 = f2->GetChisquare();
@@ -555,9 +567,16 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
                  graphErr=0;
                  delete f2;
 
-                 prChi2P2 = TMath::Prob(chi2P2, nbpf2-1);
+                 prChi2P2 = TMath::Prob(chi2P2, nbpf2Dynamic-1);
                  a2 = par[0];
                }
+             else 
+               { 
+                 unparabolicfit++;
+                 (*fFilcout) << " Warning : BP = " << busPatchId << " Manu = " << manuId <<  " Channel = " << channelId <<": parabolic fit not possible (nbpf2=" <<  nbpf2Dynamic  << ") => a2=0 and linear fit OK" << std::endl;
+                 if(unparabolicfit==1) std::cout << " Warning : BP = " << busPatchId << " Manu = " << manuId <<  ": no parabolic fit for some channels (nbpf2=" <<  nbpf2Dynamic  << "), linear fit is OK (see .log for details)" << std::endl;
+                 a2=0. ; prChi2P2=0. ;
+               }
 
              par[0] = a0;
              par[1] = a1;
@@ -600,7 +619,7 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
              q=0;  
              par[1]=0.5; a1=0.5; p1=0;
              par[2]=0.;  a2=0.;  p2=0;
-             threshold=kADCMax;        
+             threshold=ADCMax();       
 
              // bad calibration counter
              char bpmanuname[256];
@@ -629,13 +648,14 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
                {
                  //                  if(q==0  and  nplot < 100)
                  //      if(p1>1 && p2==0  and  nplot < 100)
-                           if(p1>10 && p2>10  and  nplot < 100)
+                 //                        if(p1>10 && p2>10  and  nplot < 100)
                  //    if(p1>=1 and p1<=2  and  nplot < 100)
 //               if((p1==1 || p2==1) and  nplot < 100)
+                           if(nbpf2Dynamic<nbpf2  and  nplot < 100)
                    {
                      nplot++;
                      //              cout << " nplot = " << nplot << endl;
-                     TF1 *f2Calib = new TF1("f2Calib",funcCalib,0.,kADCMax,NFITPARAMS);
+                     TF1 *f2Calib = new TF1("f2Calib",funcCalib,0.,ADCMax(),NFITPARAMS);
 
                      graphErr = new TGraphErrors(nEntries,pedMean,injCharge,pedSigma,injChargeErr);
 
@@ -657,16 +677,13 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
                      delete f2Calib;
                    }
                }
-
-
              tg->Fill();
            }
-
-
          pfilew << WriteGainData(busPatchId,manuId,channelId,par[1],par[2],threshold,q);
        }
       nmanu++;
-      if(nmanu % 500 == 0)std::cout << " Nb manu = " << nmanu << std::endl;
+      Int_t step=500;
+      if(nmanu % step == 0)std::cout << " Nb manu = " << nmanu << std::endl;
     }
 
   //      print in logfile