]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/AliMUONGain.cxx
recover the detector tag int and not the string from the AOD particle
[u/mrichter/AliRoot.git] / MUON / AliMUONGain.cxx
index 8fdeea9ca2a3ef431bd7679ce1286ca027f6fec5..8c3a0a6de9990150020736d49121937aad32c78c 100644 (file)
@@ -34,6 +34,7 @@
 #include <Riostream.h>
 
 #include <sstream>
+#include <cstdio>
 
 #define  NFITPARAMS 4
 
@@ -83,6 +84,10 @@ namespace {
   
 }
 
+using std::endl;
+using std::cout;
+using std::istringstream;
+using std::ostringstream;
 /// \cond CLASSIMP
 ClassImp(AliMUONGain)
 /// \endcond
@@ -165,7 +170,7 @@ void AliMUONGain::MakePedStoreForGain(TString shuttleFile)
       // compute and store mean DAC values (like pedestals)
       flatFile = Form("%s.ped",fPrefixDA.Data());
       outputFile=flatFile;
-      cout << "\n" << fPrefixDA.Data() << " : Flat file  generated  : " << flatFile.Data() << "\n";
+      cout << "\n" << fPrefixLDC.Data() << " : Flat file  generated  : " << flatFile.Data() << "\n";
       if (!outputFile.IsNull())  
       {
         ofstream out(outputFile.Data());
@@ -232,7 +237,6 @@ TString AliMUONGain::WriteGainHeader(Int_t nInit, Int_t nEntries, Int_t nbpf2, I
   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<<"//" << endl; 
 
   stream<<"//   RUN     DAC   " << endl;
   stream<<"//-----------------" << endl;
@@ -260,17 +264,19 @@ TString AliMUONGain::WriteGainData(Int_t BP, Int_t Manu, Int_t ch, Double_t p1,
 //_______________________________________________________________________________
 void AliMUONGain::MakeGainStore(TString shuttleFile)
 {
+  Int_t status=0;
   /// Store gains in ASCII files
   ofstream fileout;
   ofstream filcouc;
   TString tempstring;  
   TString filename; 
+  char* detail;
 
   Double_t goodA1Min =  0.5;
   Double_t goodA1Max =  2.;
 //   Double_t goodA2Min = -0.5E-03;
 //   Double_t goodA2Max =  1.E-03;
-  Double_t goodA2Min = -0.5E-01; // changed 28/10/2009 (JLC) <=> no condition on a2
+  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);
@@ -317,27 +323,23 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   Double_t injChargeErr[11];
   for ( Int_t i=0 ; i<11 ; i++) {injCharge[i]=0.;injChargeErr[i]=1.;};
 
-  // some print
-  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) {
-    cout<< " Run = " << run[i]->GetFirst() << "    DAC = " << run[i]->GetSecond() << endl;
-    numrun[i] = run[i]->GetFirst();
-    injCharge[i] = run[i]->GetSecond();
-    injChargeErr[i] = 0.01*injCharge[i];
-    if(injChargeErr[i] <= 1.) injChargeErr[i]=1.;
-  }
-  cout << "" << endl;
-
   //  print out in .log file
-
+  detail=Form("\n%s : ------  MUONTRKGAINda for Gain computing (Last Run = %d) ------\n",fPrefixLDC.Data(),fRunNumber); printf("%s",detail);
   (*fFilcout)<<"\n\n//=================================================" << endl;
   (*fFilcout)<<"//    MUONTRKGAINda: Gain Computing  Run = " << fRunNumber << endl;
   (*fFilcout)<<"//    RootDataFile  = "<< fRootDataFileName.Data() << endl;
   (*fFilcout)<<"//=================================================" << endl;
   (*fFilcout)<<"//* Date          : " << fDate->AsString("l") << "\n" << endl;
 
+  (*fFilcout) << " Entries = " << nEntries << " DAC values \n" << endl; 
+  for (Int_t i = 0; i < nEntries; ++i) {
+    (*fFilcout) << " Run = " << run[i]->GetFirst() << "    DAC = " << run[i]->GetSecond() << endl;
+    numrun[i] = run[i]->GetFirst();
+    injCharge[i] = run[i]->GetSecond();
+    injChargeErr[i] = 0.01*injCharge[i];
+    if(injChargeErr[i] <= 1.) injChargeErr[i]=1.;
+  }
+  detail=Form("%s : .... Fitting .... \n",fPrefixLDC.Data()); printf("%s",detail);
 
 
   // why 2 files ? (Ch. F.)  => second file contains detailed results
@@ -345,7 +347,8 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
     if(fPrintLevel>1)
       {
         filename=Form("%s.param",fPrefixDA.Data());
-        cout << " Second fit parameter file        = " << filename.Data() << "\n";
+       detail=Form("%s : Second fit parameter file        = %s\n",fPrefixLDC.Data(),filename.Data()); printf("%s",detail);
+       //       cout << " Second fit parameter file        = " << filename.Data() << "\n";
         pfilen = fopen (filename.Data(),"w");
 
         fprintf(pfilen,"//===================================================================\n");
@@ -368,7 +371,8 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
   if(fPrintLevel>1)
     {
       filename=Form("%s.peak",fPrefixDA.Data());
-      cout << " File containing Peak mean values = " << filename << "\n";
+      detail=Form("%s : File containing Peak mean values = %s\n",fPrefixLDC.Data(),filename.Data()); printf("%s",detail);
+      //      cout << " File containing Peak mean values = " << filename << "\n";
       pfilep = fopen (filename,"w");
 
       fprintf(pfilep,"//==============================================================================================================================\n");
@@ -398,12 +402,11 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
 
   //  plot out 
 
-  TFile* gainFile = 0x0;
   TTree* tg = 0x0;
   if(fPlotLevel>0)
     {
       fHistoFileName=Form("%s.root",fPrefixDA.Data());
-      gainFile = new TFile(fHistoFileName.Data(),"RECREATE","MUON Tracking gains");
+      new TFile(fHistoFileName.Data(),"RECREATE","MUON Tracking gains");
       tg = new TTree("tg","TTree avec class Manu_DiMu");
 
       tg->Branch("bp",&busPatchId, "busPatchId/I");
@@ -442,10 +445,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();
@@ -472,7 +477,6 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
            n++;
          }
 
-
          // print_peak_mean_values
          if(fPrintLevel>1)
            {
@@ -482,29 +486,29 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
              fprintf(pfilep,"                   sig= %9.3f%9.3f%9.3f%9.3f%9.3f%9.3f%9.3f%9.3f%9.3f%9.3f%9.3f \n",pedSigma[0],pedSigma[1],pedSigma[2],pedSigma[3],pedSigma[4],pedSigma[5],pedSigma[6],pedSigma[7],pedSigma[8],pedSigma[9],pedSigma[10]);
            }
 
-         // makegain 
-
-
+         // makegain
          // Fit Method:  Linear fit over gAlinbpf1 points + parabolic fit  over nbpf2  points) 
          // nInit=1 : 1st pt DAC=0 excluded
 
          // 1. - linear fit over gAlinbpf1 points
 
-         Double_t par[4] = {0.,0.5,0.,ADCMax()};
+         Double_t par[4] = {0.,0.5,0.,static_cast<Double_t>(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.){fitproceed=0; break;}
-             if(x[j]== ADCMax())
+             //              if(x[j]>= ADCMax())
+             if(x[j]>= adcLimit)
                {
                  if(j < nbs-1){fitproceed=0; break;}
-                 else  nbpf2Dynamic=nbpf2-1;
+                 else { nbpf2Dynamic=nbpf2-1; break;}
                }
              xErr[j] = pedSigma[k];
              y[j]    = injCharge[k];
@@ -541,8 +545,8 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
              a1 = par[1];
 
              // 2. - Translation : new origin (xLim, yLim) + parabolic fit over nbf2 points
-
-             if(nbpf2Dynamic > 1)
+             //checking:         if(busPatchId ==1841 && manuId==4)nbpf2Dynamic=2;
+             if(nbpf2Dynamic > 2)
                {
                  for (Int_t j = 0; j < nbpf2Dynamic; j++)
                    {
@@ -568,6 +572,13 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
                  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;
@@ -616,7 +627,7 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
              char bpmanuname[256];
              AliMUONErrorCounter* uncalcounter;
 
-             sprintf(bpmanuname,"bp%dmanu%d",busPatchId,manuId);
+             snprintf(bpmanuname,256,"bp%dmanu%d",busPatchId,manuId);
              if (!(uncalcounter = (AliMUONErrorCounter*)uncalBuspatchManuTable->FindObject(bpmanuname)))
                {
                  // New buspatch_manu name
@@ -650,14 +661,14 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
 
                      graphErr = new TGraphErrors(nEntries,pedMean,injCharge,pedSigma,injChargeErr);
 
-                     sprintf(graphName,"BusPatch_%d_Manu_%d_Ch_%d",busPatchId, manuId,channelId);
+                     snprintf(graphName,256,"BusPatch_%d_Manu_%d_Ch_%d",busPatchId, manuId,channelId);
 
                      graphErr->SetTitle(graphName);
                      graphErr->SetMarkerColor(3);
                      graphErr->SetMarkerStyle(12);
                      graphErr->Write(graphName);
 
-                     sprintf(graphName,"f2_BusPatch_%d_Manu_%d_Ch_%d",busPatchId, manuId,channelId);
+                     snprintf(graphName,256,"f2_BusPatch_%d_Manu_%d_Ch_%d",busPatchId, manuId,channelId);
                      f2Calib->SetTitle(graphName);
                      f2Calib->SetLineColor(4);
                      f2Calib->SetParameters(par);
@@ -668,16 +679,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)printf("%s : Nb manu = %d\n",fPrefixLDC.Data(),nmanu);
     }
 
   //      print in logfile
@@ -700,36 +708,43 @@ void AliMUONGain::MakeGainStore(TString shuttleFile)
       (*fFilcout) << " Number of bad calibrated channel = " << uncalcountertotal << endl;
        
     }
+  if(nmanu && nGoodChannel) 
+    {
+      Double_t ratio_limit=0.25;
+      Double_t ratio=float (nBadChannel)/float (nmanu*64);
+
+      detail=Form("\n%s : Nb of channels in raw data = %d (%d Manu)",fPrefixLDC.Data(),nmanu*64,nmanu); 
+      (*fFilcout) << detail ; printf("%s",detail);
+      detail=Form("\n%s : Nb of calibrated channel   = %d  (%4.2f <a1< %4.2f and %4.2f <a2< %4.2f)",fPrefixLDC.Data(),nGoodChannel,goodA1Min,goodA1Max, goodA2Min,goodA2Max);
+      (*fFilcout) << detail ; printf("%s",detail);
+      detail=Form("\n%s : Nb of uncalibrated channel = %d  [%6.4f] (%d unfitted channels) ",fPrefixLDC.Data(),nBadChannel,ratio,noFitChannel);
+      (*fFilcout) << detail ; printf("%s",detail);
+
+      if(ratio > ratio_limit) { status=-1;
+       detail=Form("\n%s : !!!!! WARNING : Nb of uncalibrated channels very large : %6.4f > %6.4f (status= %d) ",fPrefixLDC.Data(),ratio,ratio_limit,status);
+       (*fFilcout) << detail ; printf("%s",detail); 
+      }
 
-
-  (*fFilcout) << "\n Nb of channels in raw data = " << nmanu*64 << " (" << nmanu << " Manu)" <<  endl;
-  (*fFilcout) << " Nb of calibrated channel   = " << nGoodChannel << " (" << goodA1Min << "<a1<" << goodA1Max 
-             << " and " << goodA2Min << "<a2<" << goodA2Max << ") " << endl;
-  (*fFilcout) << " Nb of uncalibrated channel = " << nBadChannel << " (" << noFitChannel << " unfitted)" << endl;
-
-  cout << "\n Nb of channels in raw data = " << nmanu*64 << " (" << nmanu << " Manu)" <<  endl;
-  cout << " Nb of calibrated channel   = " << nGoodChannel << " (" << goodA1Min << "<a1<" << goodA1Max 
-       << " and " << goodA2Min << "<a2<" << goodA2Max << ") " << endl;
-  cout << " Nb of uncalibrated channel = " << nBadChannel << " (" << noFitChannel << " unfitted)" << endl;
-
-  Double_t meanA1         = sumA1/(nGoodChannel);
-  Double_t meanProbChi2   = sumProbChi2/(nGoodChannel);
-  Double_t meanA2         = sumA2/(nGoodChannel);
-  Double_t meanProbChi2P2 = sumProbChi2P2/(nGoodChannel);
-
-  Double_t capaManu = 0.2; // pF
-  (*fFilcout) << "\n linear fit   : <a1> = " << meanA1 << "\t  <gain>  = " <<  1./(meanA1*capaManu) 
-             << " mV/fC (capa= " << capaManu << " pF)" << endl;
-  (*fFilcout) <<   "        Prob(chi2)>  = " <<  meanProbChi2 << endl;
-  (*fFilcout) << "\n parabolic fit: <a2> = " << meanA2  << endl;
-  (*fFilcout) <<   "        Prob(chi2)>  = " <<  meanProbChi2P2 << "\n" << endl;
-
-  cout << "\n  <gain>  = " <<  1./(meanA1*capaManu) 
-       << " mV/fC (capa= " << capaManu << " pF)" 
-       <<  "  Prob(chi2)>  = " <<  meanProbChi2 << endl;
-  
+      Double_t meanA1         = sumA1/(nGoodChannel);
+      Double_t meanProbChi2   = sumProbChi2/(nGoodChannel);
+      Double_t meanA2         = sumA2/(nGoodChannel);
+      Double_t meanProbChi2P2 = sumProbChi2P2/(nGoodChannel);
+      Double_t capaManu = 0.2; // pF
+      (*fFilcout) << "\n linear fit   : <a1> = " << meanA1 << "     Prob(chi2)>  = " <<  meanProbChi2 <<  "    <gain>  = " <<  1./(meanA1*capaManu) 
+                 << " mV/fC (capa= " << capaManu << " pF)" << endl;
+      (*fFilcout)   <<" parabolic fit: <a2> = " << meanA2  << "    Prob(chi2)>  = " <<  meanProbChi2P2 << "\n" <<  endl;
+      detail=Form("\n%s : <gain>  = %7.5f  mV/fC (capa= %3.1f pF)  Prob(chi2) = %5.3f\n" ,fPrefixLDC.Data(),1./(meanA1*capaManu),capaManu,meanProbChi2);
+      printf("%s",detail);
+    }
+  else 
+    { status=-1;
+      detail=Form("\n%s : !!!!! ERROR :  Nb of Manu = %d or Nb calibrated channel = %d !!!!! (status= %d)\n",fPrefixLDC.Data(),nmanu,nGoodChannel,status);
+      (*fFilcout) << detail ; printf("%s",detail);
+    }
   pfilew.close();
 
   if(fPlotLevel>0){tg->Write();histoFile->Close();}
   if(fPrintLevel>1){fclose(pfilep); fclose(pfilen);}
+  SetStatusDA(status);
 }