[u/mrichter/AliRoot.git] / MUON / MUONTRKGAINda.cxx

/*
 Contact: Jean-Luc Charvet <jean-luc.charvet@cea.fr>
 Link: http://aliceinfo.cern.ch/static/Offline/dimuon/muon_html/README_Mchda
 Run Type: CALIBRATION
 DA Type: LDC
 Number of events needed: 400 events for each calibration run
 Input Files: Rawdata file (DATE format)
 Output Files: local dir (not persistent) -> MUONTRKGAINda_<run#>.par
 FXS -> run<#>_MCH_<ldc>_GAINS
 Trigger types used:
*/

/**************************************************************************
* 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.                  *
**************************************************************************/

/* $Id$ */

/*
 -------------------------------------------------------------------------
 2009-05-18 New version: MUONTRKGAINda.cxx,v 1.0
 -------------------------------------------------------------------------

 Version for MUONTRKGAINda MUON tracking
 (A. Baldisseri, J.-L. Charvet)


 Rem:  AliMUON2DMap stores all channels, even those which are not connected
 if pedMean == -1, channel not connected to a pad  


*/
extern "C" {
#include <daqDA.h>
}

#include "event.h"
#include "monitor.h"

#include <Riostream.h>
#include <stdio.h>
#include <stdlib.h>
#include <sstream>
#include <math.h> 

//AliRoot
#include "AliMUONRawStreamTrackerHP.h"
#include "AliRawReader.h"
#include "AliMUONVStore.h"
#include "AliMUON2DMap.h"
#include "AliMUONCalibParamND.h"
#include "AliMpIntPair.h"
#include "AliMpConstants.h"
#include "AliRawDataErrorLog.h"
#include "AliMUONTrackerIO.h"

//ROOT
#include "TFile.h"
#include "TSystem.h"
#include "TTree.h"
#include "TH1F.h"
#include "TString.h"
#include "TStopwatch.h"
#include "TMath.h"
#include "TTimeStamp.h"
#include "TGraphErrors.h"
#include "TF1.h"
#include "TROOT.h"
#include "TPluginManager.h"
#include "TFitter.h"
#include "TObjString.h"
#include "THashTable.h"
#include <THashList.h>
//
//AMORE
//
#ifdef ALI_AMORE
#include <AmoreDA.h>
#endif

#include "AliMUONPedestal.h"
#include "AliMUONGain.h"
#include "AliMUONErrorCounter.h"

// main routine
int main(Int_t argc, Char_t **argv) 
{

 Int_t status=0;
 TStopwatch timers;
 timers.Start(kTRUE); 

 // needed for streamer application
 gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
					"*",
					"TStreamerInfo",
					"RIO",
					"TStreamerInfo()"); 

 // needed for Minuit plugin
 gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer",
					"Minuit",
					"TMinuitMinimizer",
					"Minuit",
					"TMinuitMinimizer(const char*)");

 Char_t prefixDA[256]="MUONTRKGAINda"; // program prefix 

 Int_t skipEvents = 0;
 Int_t maxEvents  = 1000000;
 Int_t maxDateEvents  = 1000000;
 Char_t inputFile[256]="";

 Int_t  nDateEvents = 0;
 Int_t nGlitchErrors= 0;
 Int_t nParityErrors= 0;
 Int_t nPaddingErrors= 0;
 Int_t recoverParityErrors = 1;

 TString logOutputFile;

 Char_t flatFile[256]="";
 TString shuttleFile;

 Int_t nEventsRecovered = 0;
 Int_t nEvents = 0;
 UInt_t runNumber   = 0;
 Int_t  nChannel    = 0;
 ofstream filcout;
 Int_t nIndex = -1; 

 // For DA Gain
 Int_t printLevel  = 0;  // printout (up to 2)
 Int_t plotLevel  = 2;  // plotout (up to 2)
 Int_t injCharge; 

 Int_t nEntries = daqDA_ECS_getTotalIteration(); // usually = 11 = Nb of calibration runs
 Int_t nInit=1;  // = 0 all DAC values ; = 1 DAC=0 excluded (default=1)
 Int_t nbpf1=6;  // nb of points for linear fit (default=6) 


 // decode the input line
 for (Int_t i = 1; i < argc; i++) // argument 0 is the executable name
   {
     Char_t* arg;

     arg = argv[i];
     if (arg[0] != '-') 
	{
	  // If only one argument and no "-" => DA calling from ECS
	  if (argc == 2)
	    {
	      sprintf(inputFile,argv[i]);
	    }
	  continue;
	}
     switch (arg[1])
	{
	case 'f' : 
	  i++;
	  sprintf(inputFile,argv[i]);
	  break;
	case 'a' : 
	  i++;
	  shuttleFile = argv[i];
	  break;
	case 'd' :
	  i++; 
	  printLevel=atoi(argv[i]);
	  break;
	case 'g' :
	  i++; 
	  plotLevel=atoi(argv[i]);
	  break;
	case 'i' :
	  i++; 
	  nbpf1=atoi(argv[i]);
	  break;
	case 'j' :
	  i++; 
	  nInit=atoi(argv[i]);
	  break;
	case 's' :
	  i++; 
	  skipEvents=atoi(argv[i]);
	  break;
	case 'm' :
	  i++; 
	  sscanf(argv[i],"%d",&maxDateEvents);
	  break;
	case 'n' :
	  i++; 
	  sscanf(argv[i],"%d",&maxEvents);
	  break;
	case 'p' : 
	  i++;
	  sscanf(argv[i],"%d",&recoverParityErrors);
	  break;
	case 'h' :
	  i++;
	  printf("\n******************* %s usage **********************",argv[0]);
	  printf("\nOnline (called from ECS) : %s <raw data file> (no inline options)\n",argv[0]);
	  printf("\n%s -options, the available options are :",argv[0]);
	  printf("\n-h help                    (this screen)");
	  printf("\n");
	  printf("\n Input");
	  printf("\n-f <raw data file>         (default = %s)",inputFile); 
	  printf("\n");
	  printf("\n Output");
	  printf("\n-a <Flat ASCII file>       (default = %s)",shuttleFile.Data()); 
	  printf("\n");
	  printf("\n Options");
	  printf("\n-d <print level>           (default = %d)",printLevel);
	  printf("\n-g <plot level>            (default = %d)",plotLevel);
	  printf("\n-i <nb linear points>      (default = %d)",nbpf1);
	  printf("\n-j start point of fit      (default = %d)",nInit);
	  printf("\n-m <max date events>       (default = %d)",maxDateEvents);
	  printf("\n-s <skip events>           (default = %d)",skipEvents);
	  printf("\n-n <max events>            (default = %d)",maxEvents);
	  printf("\n-p <Recover parity errors> (default = %d)",recoverParityErrors);

	  printf("\n\n");
	  exit(-1);
	default :
	  printf("%s : bad argument %s (please check %s -h)\n",argv[0],argv[i],argv[0]);
	  argc = 2; exit(-1); // exit if error
	} // end of switch  
   } // end of for i  


 //Gain object
 AliMUONGain* muonGain = new AliMUONGain();
 muonGain->SetprefixDA(prefixDA);
 muonGain->SetAliRootDataFileName(); // MUONTRKGAINda_data.root  
 muonGain->SetAliPrintLevel(printLevel);
 muonGain->SetAliPlotLevel(plotLevel);

 // Reading current iteration
 nIndex = daqDA_ECS_getCurrentIteration();
 if(nIndex<0)nIndex=0; // compute gain directly from existing root data file


 if(nIndex==0) // compute gain from existing root data file: MUONTRKGAINda_data.root
   {
     sprintf(flatFile,"%s_data.log",prefixDA);
     logOutputFile=flatFile;
     filcout.open(logOutputFile.Data());
     muonGain->SetAlifilcout(&filcout);

     sprintf(flatFile,"%s_data.par",prefixDA);
     if(shuttleFile.IsNull())shuttleFile=flatFile;

     muonGain->SetAliIndex(nIndex); // fIndex 
     muonGain->SetAliInit(nInit); // fnInit
     muonGain->SetAliNbpf1(nbpf1); // fnbpf1
     muonGain->MakeGainStore(shuttleFile);
   }


 if(nIndex>0) // normal case: reading calibration runs before computing gains
   {
     UShort_t manuId;  
     UChar_t channelId;
     UShort_t charge;
	  
     // DAC values stored in array vDAC (reading dbfile in DETDB)
     //   Int_t vDAC[11] = {0,200,400,800,1200,1600,2000,2500,3000,3500,4000}; // DAC values
     Int_t vDAC[11]; // DAC values
     Char_t *dbfile;
     dbfile="mutrkcalibvalues";
     cout << "\n *** Copy: " << dbfile << " from DetDB to working directory  *** \n" << endl;
     status=daqDA_DB_getFile(dbfile,dbfile);
     ifstream filein(dbfile,ios::in); Int_t k=0, kk;
     while (k<nEntries ) { filein >> kk >> vDAC[k] ; k++; }
		
     injCharge=vDAC[nIndex-1];

     // Rawdeader, RawStreamHP
     AliRawReader* rawReader = AliRawReader::Create(inputFile);
     AliMUONRawStreamTrackerHP* rawStream  = new AliMUONRawStreamTrackerHP(rawReader);    
     rawStream->DisableWarnings();
     rawStream->EnabbleErrorLogger();

     cout << "\n" << prefixDA << " : Reading data from file " << inputFile  << endl;

     while (rawReader->NextEvent())
	{
	  if (nDateEvents >= maxDateEvents) break;
	  if (nEvents >= maxEvents) break;
	  if (nDateEvents>0 &&  nDateEvents % 100 == 0) 	
	    cout<<"Cumulated:  DATE events = " << nDateEvents << "   Used events = " << nEvents << endl;

	  // check shutdown condition 
	  if (daqDA_checkShutdown()) 
	    break;

	  //Skip events
	  while (skipEvents)
	    {
	      rawReader->NextEvent();
	      skipEvents--;
	    }

	  Int_t eventType = rawReader->GetType();
	  runNumber = rawReader->GetRunNumber();

	  // Output log file initialisations
	  if(nDateEvents==0)
	    {
	      sprintf(flatFile,"%s_%d.log",prefixDA,runNumber);
	      logOutputFile=flatFile;

	      filcout.open(logOutputFile.Data());
	      filcout<<"//=================================================" << endl;
	      filcout<<"//       " << prefixDA << " for run = " << runNumber << "  (DAC=" << injCharge << ")" << endl;
	      filcout<<"//=================================================" << endl;
	      filcout<<"//   * Date          : " << muonGain->GetDate()->AsString("l") << "\n" << endl;

	      cout<<"\n ********  " << prefixDA << " for run = " << runNumber << "  (Index= " << nIndex << "/" << nEntries << "  DAC=" << injCharge << ") ********\n" << endl;
	      cout<<" * Date : " << muonGain->GetDate()->AsString("l") << "\n" << endl;
	    }

	  muonGain->SetAlifilcout(&filcout);

	  nDateEvents++;
	  if (eventType != PHYSICS_EVENT)
	    continue; // for the moment

	  // First lopp over DDL's to find good events
	  // Error counters per event (counters in the decoding lib are for each DDL)
	  Bool_t eventIsErrorMessage = kFALSE;
	  int eventGlitchErrors = 0;
	  int eventParityErrors = 0;
	  int eventPaddingErrors = 0;
	  rawStream->First();
	  do
	    {
	      if (rawStream->IsErrorMessage()) eventIsErrorMessage = kTRUE;
	      eventGlitchErrors += rawStream->GetGlitchErrors();
	      eventParityErrors += rawStream->GetParityErrors();
	      eventPaddingErrors += rawStream->GetPaddingErrors();
	    } while(rawStream->NextDDL()); 

	  AliMUONRawStreamTrackerHP::AliBusPatch* busPatch;
	  if (!eventIsErrorMessage) 
	    {
	      // Good events (no error) -> compute pedestal for all channels
	      rawStream->First(); 
	      while( (busPatch = (AliMUONRawStreamTrackerHP::AliBusPatch*) rawStream->Next())) 
		{
		  for(int i = 0; i < busPatch->GetLength(); ++i)
		    {
		      if (nEvents == 0) nChannel++;
		      busPatch->GetData(i, manuId, channelId, charge);
		      muonGain->MakePed(busPatch->GetBusPatchId(), (Int_t)manuId, (Int_t)channelId, (Int_t)charge);
		    }
		}
	      nEvents++;
	    }
	  else
	    {
	      // Events with errors
	      if (recoverParityErrors && eventParityErrors && !eventGlitchErrors&& !eventPaddingErrors)
		{
		  // Recover parity errors -> compute pedestal for all good buspatches
		  if ( TEST_SYSTEM_ATTRIBUTE( rawReader->GetAttributes(),
					      ATTR_ORBIT_BC )) 
		    {
		      filcout <<"Event recovered -> Period:"<<EVENT_ID_GET_PERIOD( rawReader->GetEventId() )
				  <<" Orbit:"<<EVENT_ID_GET_ORBIT( rawReader->GetEventId() )
				  <<" BunchCrossing:"<<EVENT_ID_GET_BUNCH_CROSSING( rawReader->GetEventId() )<<endl;				
		    } 
		  else 
		    {
		      filcout <<"Event recovered -> nbInRun:"<<EVENT_ID_GET_NB_IN_RUN( rawReader->GetEventId() )
				  <<" burstNb:"<<EVENT_ID_GET_BURST_NB( rawReader->GetEventId() )
				  <<" nbInBurst:"<<EVENT_ID_GET_NB_IN_BURST( rawReader->GetEventId() )<<endl;
		    }
		  rawStream->First();
		  while( (busPatch = (AliMUONRawStreamTrackerHP::AliBusPatch*) rawStream->Next())) 
		    {
		      // Check the buspatch -> if error not use it in the pedestal calculation
		      int errorCount = 0;
		      for(int i = 0; i < busPatch->GetLength(); ++i)
			{
			  if (!busPatch->IsParityOk(i)) errorCount++;
			}
		      if (!errorCount) 
			{
			  // Good buspatch
			  for(int i = 0; i < busPatch->GetLength(); ++i)
			    {
			      if (nEvents == 0) nChannel++;
			      busPatch->GetData(i, manuId, channelId, charge);
			      // if (busPatch->GetBusPatchId()==1719 && manuId == 1 && channelId == 0) cout <<"Recovered charge "<<charge<<endl;
			      muonGain->MakePed(busPatch->GetBusPatchId(), (Int_t)manuId, (Int_t)channelId, (Int_t)charge);
			    }
			}
		      else
			{
			  char bpname[256];
			  AliMUONErrorCounter* errorCounter;
			  // Bad buspatch -> not used (just print)
			  filcout<<"bpId "<<busPatch->GetBusPatchId()<<" words "<<busPatch->GetLength()
				     <<" parity errors "<<errorCount<<endl;
			  // Number of events where this buspatch is missing
			  sprintf(bpname,"bp%d",busPatch->GetBusPatchId());						
			  if (!(errorCounter = (AliMUONErrorCounter*) (muonGain->GetErrorBuspatchTable()->FindObject(bpname))))
			    {
			      // New buspatch
			      errorCounter = new AliMUONErrorCounter(busPatch->GetBusPatchId());
			      errorCounter->SetName(bpname);
			      muonGain->GetErrorBuspatchTable()->Add(errorCounter);
			    }
			  else
			    {
			      // Existing buspatch
			      errorCounter->Increment();
			    }	
			  // errorCounter->Print();						
			} // end of if (!errorCount)
		    } // end of while( (busPatch = (AliMUONRawStreamTrackerHP ...
		  nEvents++;
		  nEventsRecovered++;
		} //end of if (recoverParityErrors && eventParityErrors && !eventGlitchErrors&& !eventPaddingErrors)
	      else
		{
		  // Fatal errors reject the event
		  if ( TEST_SYSTEM_ATTRIBUTE( rawReader->GetAttributes(),
					      ATTR_ORBIT_BC )) 
		    {
		      filcout <<"Event rejected -> Period:"<<EVENT_ID_GET_PERIOD( rawReader->GetEventId() )
				  <<" Orbit:"<<EVENT_ID_GET_ORBIT( rawReader->GetEventId() )
				  <<" BunchCrossing:"<<EVENT_ID_GET_BUNCH_CROSSING( rawReader->GetEventId() )<<endl;				
		    } 
		  else 
		    {
		      filcout <<"Event rejected -> nbInRun:"<<EVENT_ID_GET_NB_IN_RUN( rawReader->GetEventId() )
				  <<" burstNb:"<<EVENT_ID_GET_BURST_NB( rawReader->GetEventId() )
				  <<" nbInBurst:"<<EVENT_ID_GET_NB_IN_BURST( rawReader->GetEventId() )<<endl;

		    }
		} // end of if (!rawStream->GetGlitchErrors() && !rawStream->GetPaddingErrors() ...
	      filcout<<"Number of errors : Glitch "<<eventGlitchErrors
			 <<" Parity "<<eventParityErrors
			 <<" Padding "<<eventPaddingErrors<<endl;
	      filcout<<endl;			
	    } // end of if (!rawStream->IsErrorMessage())

	  if (eventGlitchErrors)  nGlitchErrors++;
	  if (eventParityErrors)  nParityErrors++;
	  if (eventPaddingErrors) nPaddingErrors++;

	} // while (rawReader->NextEvent())
     delete rawReader;
     delete rawStream;

     // process and store mean peak values in .root file
    sprintf(flatFile,"%s_%d.par",prefixDA,runNumber);
     if(shuttleFile.IsNull())shuttleFile=flatFile;
     injCharge=vDAC[nIndex-1];
     muonGain->SetAliIndex(nIndex); // fIndex 
     muonGain->SetAliInjCharge(injCharge);
     muonGain->SetAliNEvents(nEvents);
     muonGain->SetAliRunNumber(runNumber);
     muonGain->SetAliNChannel(nChannel);
     muonGain->MakePedStoreForGain(shuttleFile);

     // writing some counters
     cout << endl;
     cout << prefixDA << " : Nb of DATE events           = " << nDateEvents    << endl;
     cout << prefixDA << " : Nb of Glitch errors         = "   << nGlitchErrors  << endl;
     cout << prefixDA << " : Nb of Parity errors         = "   << nParityErrors  << endl;
     cout << prefixDA << " : Nb of Padding errors        = "   << nPaddingErrors << endl;		
     cout << prefixDA << " : Nb of events recovered      = "   << nEventsRecovered<< endl;
     cout << prefixDA << " : Nb of events without errors = "   << nEvents-nEventsRecovered<< endl;
     cout << prefixDA << " : Nb of events used           = "   << nEvents        << endl;

     filcout << endl;
     filcout << prefixDA << " : Nb of DATE events           = " << nDateEvents    << endl;
     filcout << prefixDA << " : Nb of Glitch errors         = "   << nGlitchErrors << endl;
     filcout << prefixDA << " : Nb of Parity errors         = "   << nParityErrors << endl;
     filcout << prefixDA << " : Nb of Padding errors        = "   << nPaddingErrors << endl;
     filcout << prefixDA << " : Nb of events recovered      = "   << nEventsRecovered<< endl;	
     filcout << prefixDA << " : Nb of events without errors = "   << nEvents-nEventsRecovered<< endl;
     filcout << prefixDA << " : Nb of events used           = "   << nEvents        << endl;

     // Computing gain 
     if(nIndex==nEntries)
	{
	  muonGain->SetAliInit(nInit); // fnInit
	  muonGain->SetAliEntries(nEntries); // fnEntries
	  muonGain->SetAliNbpf1(nbpf1); // fnbpf1
     	  muonGain->MakeGainStore(shuttleFile);
	}

//       Copying files to local DB folder defined by DAQ_DETDB_LOCAL
	Char_t *dir;
      dir= getenv("DAQ_DETDB_LOCAL");
      unsigned int nLastVersions = 2;
      cout << "\n *** Output files stored locally in " << dir << " (nb of previous versions = " << nLastVersions << ") ***" << endl;
      filcout << "\n *** Output files stored locally in " << dir << " (nb of previous versions = " << nLastVersions << ") ***" << endl;

      filcout << endl;
      filcout << prefixDA << " : Root data file         : " << muonGain->GetRootDataFileName() << endl;
      filcout << prefixDA << " : Output logfile         : " << logOutputFile  << endl;
      filcout << prefixDA << " : Gain Histo file        : " << muonGain->GetHistoFileName() << endl;
      filcout << prefixDA << " : Gain file (to SHUTTLE) : " << shuttleFile << endl;

      status = daqDA_localDB_storeFile(logOutputFile.Data(),nLastVersions);
      if(nIndex==nEntries)status = daqDA_localDB_storeFile(muonGain->GetRootDataFileName(),nLastVersions);
      status = daqDA_localDB_storeFile(muonGain->GetHistoFileName(),nLastVersions);
      status = daqDA_localDB_storeFile(shuttleFile.Data(),nLastVersions);

    }

// ouput files
  cout << endl;
  cout << prefixDA << " : Root data file         : " << muonGain->GetRootDataFileName() << endl;
  cout << prefixDA << " : Output logfile         : " << logOutputFile  << endl;
  cout << prefixDA << " : Gain Histo file        : " << muonGain->GetHistoFileName() << endl;
  cout << prefixDA << " : Gain file (to SHUTTLE) : " << shuttleFile << endl;   
 

  // Transferring to OCDB via the SHUTTLE
  printf("\n *****  STORE FILE in FES ****** \n");

  // be sure that env variable DAQDALIB_PATH is set in script file
  //       gSystem->Setenv("DAQDALIB_PATH", "$DATE_SITE/infoLogger");

  status = daqDA_FES_storeFile(shuttleFile.Data(),"GAINS");
  if (status) 
    {
      printf(" Failed to export file : %d\n",status);
    }
  else printf(" %s successfully exported to FES  \n",shuttleFile.Data());


  filcout.close();
  timers.Stop();
  printf("\nExecution time : R:%7.2fs C:%7.2fs\n", timers.RealTime(), timers.CpuTime());
  return status;
}
Commit	Line	Data
dc928fbf	1	/*
1dd09aa3	2	Contact: Jean-Luc Charvet <jean-luc.charvet@cea.fr>
	3	Link: http://aliceinfo.cern.ch/static/Offline/dimuon/muon_html/README_Mchda
	4	Run Type: CALIBRATION
	5	DA Type: LDC
	6	Number of events needed: 400 events for each calibration run
	7	Input Files: Rawdata file (DATE format)
	8	Output Files: local dir (not persistent) -> MUONTRKGAINda_<run#>.par
	9	FXS -> run<#>_MCH_<ldc>_GAINS
	10	Trigger types used:
dc928fbf	11	*/
	12
	13	/**************************************************************************
1dd09aa3	14	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	15	* *
	16	* Author: The ALICE Off-line Project. *
	17	* Contributors are mentioned in the code where appropriate. *
	18	* *
	19	* Permission to use, copy, modify and distribute this software and its *
	20	* documentation strictly for non-commercial purposes is hereby granted *
	21	* without fee, provided that the above copyright notice appears in all *
	22	* copies and that both the copyright notice and this permission notice *
	23	* appear in the supporting documentation. The authors make no claims *
	24	* about the suitability of this software for any purpose. It is *
	25	* provided "as is" without express or implied warranty. *
	26	**************************************************************************/
dc928fbf	27
	28	/* $Id$ */
	29
	30	/*
1dd09aa3	31	-------------------------------------------------------------------------
	32	2009-05-18 New version: MUONTRKGAINda.cxx,v 1.0
	33	-------------------------------------------------------------------------
dc928fbf	34
1dd09aa3	35	Version for MUONTRKGAINda MUON tracking
1dd09aa3	36	(A. Baldisseri, J.-L. Charvet)
dc928fbf	37
dc928fbf	38
1dd09aa3	39	Rem: AliMUON2DMap stores all channels, even those which are not connected
1dd09aa3	40	if pedMean == -1, channel not connected to a pad
dc928fbf	41
	42
	43	*/
	44	extern "C" {
	45	#include <daqDA.h>
	46	}
	47
	48	#include "event.h"
	49	#include "monitor.h"
	50
	51	#include <Riostream.h>
	52	#include <stdio.h>
	53	#include <stdlib.h>
	54	#include <sstream>
	55	#include <math.h>
	56
	57	//AliRoot
	58	#include "AliMUONRawStreamTrackerHP.h"
	59	#include "AliRawReader.h"
	60	#include "AliMUONVStore.h"
	61	#include "AliMUON2DMap.h"
	62	#include "AliMUONCalibParamND.h"
	63	#include "AliMpIntPair.h"
	64	#include "AliMpConstants.h"
	65	#include "AliRawDataErrorLog.h"
	66	#include "AliMUONTrackerIO.h"
	67
	68	//ROOT
	69	#include "TFile.h"
	70	#include "TSystem.h"
	71	#include "TTree.h"
	72	#include "TH1F.h"
	73	#include "TString.h"
	74	#include "TStopwatch.h"
	75	#include "TMath.h"
	76	#include "TTimeStamp.h"
	77	#include "TGraphErrors.h"
	78	#include "TF1.h"
	79	#include "TROOT.h"
	80	#include "TPluginManager.h"
	81	#include "TFitter.h"
	82	#include "TObjString.h"
	83	#include "THashTable.h"
	84	#include <THashList.h>
	85	//
	86	//AMORE
	87	//
	88	#ifdef ALI_AMORE
	89	#include <AmoreDA.h>
	90	#endif
	91
	92	#include "AliMUONPedestal.h"
	93	#include "AliMUONGain.h"
	94	#include "AliMUONErrorCounter.h"
	95
	96	// main routine
	97	int main(Int_t argc, Char_t **argv)
	98	{
	99
1dd09aa3	100	Int_t status=0;
	101	TStopwatch timers;
	102	timers.Start(kTRUE);
dc928fbf	103
1dd09aa3	104	// needed for streamer application
1dd09aa3	105	gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
dc928fbf	106	"*",
	107	"TStreamerInfo",
	108	"RIO",
	109	"TStreamerInfo()");
	110
1dd09aa3	111	// needed for Minuit plugin
1dd09aa3	112	gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer",
dc928fbf	113	"Minuit",
	114	"TMinuitMinimizer",
	115	"Minuit",
	116	"TMinuitMinimizer(const char*)");
	117
1dd09aa3	118	Char_t prefixDA[256]="MUONTRKGAINda"; // program prefix
dc928fbf	119
1dd09aa3	120	Int_t skipEvents = 0;
	121	Int_t maxEvents = 1000000;
	122	Int_t maxDateEvents = 1000000;
	123	Char_t inputFile[256]="";
dc928fbf	124
1dd09aa3	125	Int_t nDateEvents = 0;
	126	Int_t nGlitchErrors= 0;
	127	Int_t nParityErrors= 0;
	128	Int_t nPaddingErrors= 0;
	129	Int_t recoverParityErrors = 1;
dc928fbf	130
1dd09aa3	131	TString logOutputFile;
dc928fbf	132
1dd09aa3	133	Char_t flatFile[256]="";
1dd09aa3	134	TString shuttleFile;
dc928fbf	135
1dd09aa3	136	Int_t nEventsRecovered = 0;
	137	Int_t nEvents = 0;
	138	UInt_t runNumber = 0;
	139	Int_t nChannel = 0;
	140	ofstream filcout;
	141	Int_t nIndex = -1;
dc928fbf	142
1dd09aa3	143	// For DA Gain
	144	Int_t printLevel = 0; // printout (up to 2)
	145	Int_t plotLevel = 2; // plotout (up to 2)
	146	Int_t injCharge;
dc928fbf	147
1dd09aa3	148	Int_t nEntries = daqDA_ECS_getTotalIteration(); // usually = 11 = Nb of calibration runs
	149	Int_t nInit=1; // = 0 all DAC values ; = 1 DAC=0 excluded (default=1)
	150	Int_t nbpf1=6; // nb of points for linear fit (default=6)
dc928fbf	151
dc928fbf	152
1dd09aa3	153	// decode the input line
	154	for (Int_t i = 1; i < argc; i++) // argument 0 is the executable name
	155	{
	156	Char_t* arg;
dc928fbf	157
1dd09aa3	158	arg = argv[i];
1dd09aa3	159	if (arg[0] != '-')
dc928fbf	160	{
	161	// If only one argument and no "-" => DA calling from ECS
	162	if (argc == 2)
	163	{
	164	sprintf(inputFile,argv[i]);
	165	}
	166	continue;
	167	}
1dd09aa3	168	switch (arg[1])
dc928fbf	169	{
	170	case 'f' :
	171	i++;
	172	sprintf(inputFile,argv[i]);
	173	break;
	174	case 'a' :
	175	i++;
	176	shuttleFile = argv[i];
	177	break;
	178	case 'd' :
	179	i++;
	180	printLevel=atoi(argv[i]);
	181	break;
	182	case 'g' :
	183	i++;
	184	plotLevel=atoi(argv[i]);
	185	break;
	186	case 'i' :
	187	i++;
	188	nbpf1=atoi(argv[i]);
	189	break;
	190	case 'j' :
	191	i++;
	192	nInit=atoi(argv[i]);
	193	break;
	194	case 's' :
	195	i++;
	196	skipEvents=atoi(argv[i]);
	197	break;
	198	case 'm' :
	199	i++;
	200	sscanf(argv[i],"%d",&maxDateEvents);
	201	break;
	202	case 'n' :
	203	i++;
	204	sscanf(argv[i],"%d",&maxEvents);
	205	break;
	206	case 'p' :
	207	i++;
	208	sscanf(argv[i],"%d",&recoverParityErrors);
	209	break;
	210	case 'h' :
	211	i++;
	212	printf("\n***************** %s usage ********************",argv[0]);
	213	printf("\nOnline (called from ECS) : %s <raw data file> (no inline options)\n",argv[0]);
	214	printf("\n%s -options, the available options are :",argv[0]);
	215	printf("\n-h help (this screen)");
	216	printf("\n");
	217	printf("\n Input");
	218	printf("\n-f <raw data file> (default = %s)",inputFile);
	219	printf("\n");
	220	printf("\n Output");
	221	printf("\n-a <Flat ASCII file> (default = %s)",shuttleFile.Data());
	222	printf("\n");
	223	printf("\n Options");
	224	printf("\n-d <print level> (default = %d)",printLevel);
	225	printf("\n-g <plot level> (default = %d)",plotLevel);
	226	printf("\n-i <nb linear points> (default = %d)",nbpf1);
	227	printf("\n-j start point of fit (default = %d)",nInit);
	228	printf("\n-m <max date events> (default = %d)",maxDateEvents);
	229	printf("\n-s <skip events> (default = %d)",skipEvents);
	230	printf("\n-n <max events> (default = %d)",maxEvents);
	231	printf("\n-p <Recover parity errors> (default = %d)",recoverParityErrors);
	232
233	printf("\n\n");
234	exit(-1);
235	default :
236	printf("%s : bad argument %s (please check %s -h)\n",argv[0],argv[i],argv[0]);
237	argc = 2; exit(-1); // exit if error
238	} // end of switch
1dd09aa3	239	} // end of for i
dc928fbf	240
dc928fbf	241
1dd09aa3	242	//Gain object
	243	AliMUONGain* muonGain = new AliMUONGain();
	244	muonGain->SetprefixDA(prefixDA);
	245	muonGain->SetAliRootDataFileName(); // MUONTRKGAINda_data.root
	246	muonGain->SetAliPrintLevel(printLevel);
	247	muonGain->SetAliPlotLevel(plotLevel);
dc928fbf	248
1dd09aa3	249	// Reading current iteration
	250	nIndex = daqDA_ECS_getCurrentIteration();
	251	if(nIndex<0)nIndex=0; // compute gain directly from existing root data file
dc928fbf	252
	253
	254
	255
1dd09aa3	256	if(nIndex==0) // compute gain from existing root data file: MUONTRKGAINda_data.root
	257	{
	258	sprintf(flatFile,"%s_data.log",prefixDA);
	259	logOutputFile=flatFile;
	260	filcout.open(logOutputFile.Data());
	261	muonGain->SetAlifilcout(&filcout);
dc928fbf	262
1dd09aa3	263	sprintf(flatFile,"%s_data.par",prefixDA);
1dd09aa3	264	if(shuttleFile.IsNull())shuttleFile=flatFile;
dc928fbf	265
1dd09aa3	266	muonGain->SetAliIndex(nIndex); // fIndex
	267	muonGain->SetAliInit(nInit); // fnInit
	268	muonGain->SetAliNbpf1(nbpf1); // fnbpf1
	269	muonGain->MakeGainStore(shuttleFile);
	270	}
dc928fbf	271
1dd09aa3	272
	273
	274	if(nIndex>0) // normal case: reading calibration runs before computing gains
	275	{
	276	UShort_t manuId;
	277	UChar_t channelId;
	278	UShort_t charge;
dc928fbf	279
1dd09aa3	280	// DAC values stored in array vDAC (reading dbfile in DETDB)
	281	// Int_t vDAC[11] = {0,200,400,800,1200,1600,2000,2500,3000,3500,4000}; // DAC values
	282	Int_t vDAC[11]; // DAC values
	283	Char_t *dbfile;
	284	dbfile="mutrkcalibvalues";
	285	cout << "\n * Copy: " << dbfile << " from DetDB to working directory * \n" << endl;
	286	status=daqDA_DB_getFile(dbfile,dbfile);
	287	ifstream filein(dbfile,ios::in); Int_t k=0, kk;
	288	while (k<nEntries ) { filein >> kk >> vDAC[k] ; k++; }
dc928fbf	289
1dd09aa3	290	injCharge=vDAC[nIndex-1];
dc928fbf	291
1dd09aa3	292	// Rawdeader, RawStreamHP
	293	AliRawReader* rawReader = AliRawReader::Create(inputFile);
	294	AliMUONRawStreamTrackerHP* rawStream = new AliMUONRawStreamTrackerHP(rawReader);
	295	rawStream->DisableWarnings();
	296	rawStream->EnabbleErrorLogger();
dc928fbf	297
1dd09aa3	298	cout << "\n" << prefixDA << " : Reading data from file " << inputFile << endl;
dc928fbf	299
1dd09aa3	300	while (rawReader->NextEvent())
dc928fbf	301	{
	302	if (nDateEvents >= maxDateEvents) break;
	303	if (nEvents >= maxEvents) break;
	304	if (nDateEvents>0 && nDateEvents % 100 == 0)
	305	cout<<"Cumulated: DATE events = " << nDateEvents << " Used events = " << nEvents << endl;
	306
	307	// check shutdown condition
	308	if (daqDA_checkShutdown())
	309	break;
	310
	311	//Skip events
	312	while (skipEvents)
	313	{
	314	rawReader->NextEvent();
	315	skipEvents--;
	316	}
	317
	318	Int_t eventType = rawReader->GetType();
	319	runNumber = rawReader->GetRunNumber();
	320
	321	// Output log file initialisations
	322	if(nDateEvents==0)
	323	{
	324	sprintf(flatFile,"%s_%d.log",prefixDA,runNumber);
	325	logOutputFile=flatFile;
	326
	327	filcout.open(logOutputFile.Data());
	328	filcout<<"//=================================================" << endl;
	329	filcout<<"// " << prefixDA << " for run = " << runNumber << " (DAC=" << injCharge << ")" << endl;
	330	filcout<<"//=================================================" << endl;
	331	filcout<<"// * Date : " << muonGain->GetDate()->AsString("l") << "\n" << endl;
	332
	333	cout<<"\n ****** " << prefixDA << " for run = " << runNumber << " (Index= " << nIndex << "/" << nEntries << " DAC=" << injCharge << ") ******\n" << endl;
	334	cout<<" * Date : " << muonGain->GetDate()->AsString("l") << "\n" << endl;
	335	}
	336
	337	muonGain->SetAlifilcout(&filcout);
	338
	339	nDateEvents++;
	340	if (eventType != PHYSICS_EVENT)
	341	continue; // for the moment
	342
	343	// First lopp over DDL's to find good events
	344	// Error counters per event (counters in the decoding lib are for each DDL)
	345	Bool_t eventIsErrorMessage = kFALSE;
	346	int eventGlitchErrors = 0;
	347	int eventParityErrors = 0;
	348	int eventPaddingErrors = 0;
	349	rawStream->First();
	350	do
	351	{
	352	if (rawStream->IsErrorMessage()) eventIsErrorMessage = kTRUE;
	353	eventGlitchErrors += rawStream->GetGlitchErrors();
	354	eventParityErrors += rawStream->GetParityErrors();
	355	eventPaddingErrors += rawStream->GetPaddingErrors();
	356	} while(rawStream->NextDDL());
	357
	358	AliMUONRawStreamTrackerHP::AliBusPatch* busPatch;
	359	if (!eventIsErrorMessage)
	360	{
	361	// Good events (no error) -> compute pedestal for all channels
	362	rawStream->First();
	363	while( (busPatch = (AliMUONRawStreamTrackerHP::AliBusPatch*) rawStream->Next()))
	364	{
365	for(int i = 0; i < busPatch->GetLength(); ++i)
366	{
367	if (nEvents == 0) nChannel++;
368	busPatch->GetData(i, manuId, channelId, charge);
369	muonGain->MakePed(busPatch->GetBusPatchId(), (Int_t)manuId, (Int_t)channelId, (Int_t)charge);
370	}
371	}
372	nEvents++;
373	}
374	else
375	{
376	// Events with errors
377	if (recoverParityErrors && eventParityErrors && !eventGlitchErrors&& !eventPaddingErrors)
378	{
379	// Recover parity errors -> compute pedestal for all good buspatches
380	if ( TEST_SYSTEM_ATTRIBUTE( rawReader->GetAttributes(),
381	ATTR_ORBIT_BC ))
382	{
383	filcout <<"Event recovered -> Period:"<<EVENT_ID_GET_PERIOD( rawReader->GetEventId() )
384	<<" Orbit:"<<EVENT_ID_GET_ORBIT( rawReader->GetEventId() )
385	<<" BunchCrossing:"<<EVENT_ID_GET_BUNCH_CROSSING( rawReader->GetEventId() )<<endl;
386	}
387	else
388	{
389	filcout <<"Event recovered -> nbInRun:"<<EVENT_ID_GET_NB_IN_RUN( rawReader->GetEventId() )
390	<<" burstNb:"<<EVENT_ID_GET_BURST_NB( rawReader->GetEventId() )
391	<<" nbInBurst:"<<EVENT_ID_GET_NB_IN_BURST( rawReader->GetEventId() )<<endl;
392	}
393	rawStream->First();
394	while( (busPatch = (AliMUONRawStreamTrackerHP::AliBusPatch*) rawStream->Next()))
395	{
396	// Check the buspatch -> if error not use it in the pedestal calculation
397	int errorCount = 0;
398	for(int i = 0; i < busPatch->GetLength(); ++i)
399	{
400	if (!busPatch->IsParityOk(i)) errorCount++;
401	}
402	if (!errorCount)
403	{
404	// Good buspatch
405	for(int i = 0; i < busPatch->GetLength(); ++i)
406	{
407	if (nEvents == 0) nChannel++;
408	busPatch->GetData(i, manuId, channelId, charge);
409	// if (busPatch->GetBusPatchId()==1719 && manuId == 1 && channelId == 0) cout <<"Recovered charge "<<charge<<endl;
410	muonGain->MakePed(busPatch->GetBusPatchId(), (Int_t)manuId, (Int_t)channelId, (Int_t)charge);
411	}
412	}
413	else
414	{
415	char bpname[256];
416	AliMUONErrorCounter* errorCounter;
417	// Bad buspatch -> not used (just print)
418	filcout<<"bpId "<<busPatch->GetBusPatchId()<<" words "<<busPatch->GetLength()
419	<<" parity errors "<<errorCount<<endl;
420	// Number of events where this buspatch is missing
421	sprintf(bpname,"bp%d",busPatch->GetBusPatchId());
422	if (!(errorCounter = (AliMUONErrorCounter*) (muonGain->GetErrorBuspatchTable()->FindObject(bpname))))
423	{
424	// New buspatch
425	errorCounter = new AliMUONErrorCounter(busPatch->GetBusPatchId());
426	errorCounter->SetName(bpname);
427	muonGain->GetErrorBuspatchTable()->Add(errorCounter);
428	}
429	else
430	{
431	// Existing buspatch
432	errorCounter->Increment();
433	}
434	// errorCounter->Print();
435	} // end of if (!errorCount)
436	} // end of while( (busPatch = (AliMUONRawStreamTrackerHP ...
437	nEvents++;
438	nEventsRecovered++;
439	} //end of if (recoverParityErrors && eventParityErrors && !eventGlitchErrors&& !eventPaddingErrors)
440	else
441	{
442	// Fatal errors reject the event
443	if ( TEST_SYSTEM_ATTRIBUTE( rawReader->GetAttributes(),
444	ATTR_ORBIT_BC ))
445	{
446	filcout <<"Event rejected -> Period:"<<EVENT_ID_GET_PERIOD( rawReader->GetEventId() )
447	<<" Orbit:"<<EVENT_ID_GET_ORBIT( rawReader->GetEventId() )
448	<<" BunchCrossing:"<<EVENT_ID_GET_BUNCH_CROSSING( rawReader->GetEventId() )<<endl;
449	}
450	else
451	{
452	filcout <<"Event rejected -> nbInRun:"<<EVENT_ID_GET_NB_IN_RUN( rawReader->GetEventId() )
453	<<" burstNb:"<<EVENT_ID_GET_BURST_NB( rawReader->GetEventId() )
454	<<" nbInBurst:"<<EVENT_ID_GET_NB_IN_BURST( rawReader->GetEventId() )<<endl;
455
456	}
457	} // end of if (!rawStream->GetGlitchErrors() && !rawStream->GetPaddingErrors() ...
458	filcout<<"Number of errors : Glitch "<<eventGlitchErrors
459	<<" Parity "<<eventParityErrors
460	<<" Padding "<<eventPaddingErrors<<endl;
461	filcout<<endl;
462	} // end of if (!rawStream->IsErrorMessage())
463
464	if (eventGlitchErrors) nGlitchErrors++;
465	if (eventParityErrors) nParityErrors++;
466	if (eventPaddingErrors) nPaddingErrors++;
467
468	} // while (rawReader->NextEvent())
1dd09aa3	469	delete rawReader;
	470	delete rawStream;
	471
	472	// process and store mean peak values in .root file
	473	sprintf(flatFile,"%s_%d.par",prefixDA,runNumber);
	474	if(shuttleFile.IsNull())shuttleFile=flatFile;
	475	injCharge=vDAC[nIndex-1];
	476	muonGain->SetAliIndex(nIndex); // fIndex
	477	muonGain->SetAliInjCharge(injCharge);
	478	muonGain->SetAliNEvents(nEvents);
	479	muonGain->SetAliRunNumber(runNumber);
	480	muonGain->SetAliNChannel(nChannel);
	481	muonGain->MakePedStoreForGain(shuttleFile);
	482
	483	// writing some counters
	484	cout << endl;
	485	cout << prefixDA << " : Nb of DATE events = " << nDateEvents << endl;
	486	cout << prefixDA << " : Nb of Glitch errors = " << nGlitchErrors << endl;
	487	cout << prefixDA << " : Nb of Parity errors = " << nParityErrors << endl;
	488	cout << prefixDA << " : Nb of Padding errors = " << nPaddingErrors << endl;
	489	cout << prefixDA << " : Nb of events recovered = " << nEventsRecovered<< endl;
	490	cout << prefixDA << " : Nb of events without errors = " << nEvents-nEventsRecovered<< endl;
	491	cout << prefixDA << " : Nb of events used = " << nEvents << endl;
	492
	493	filcout << endl;
	494	filcout << prefixDA << " : Nb of DATE events = " << nDateEvents << endl;
	495	filcout << prefixDA << " : Nb of Glitch errors = " << nGlitchErrors << endl;
	496	filcout << prefixDA << " : Nb of Parity errors = " << nParityErrors << endl;
	497	filcout << prefixDA << " : Nb of Padding errors = " << nPaddingErrors << endl;
	498	filcout << prefixDA << " : Nb of events recovered = " << nEventsRecovered<< endl;
	499	filcout << prefixDA << " : Nb of events without errors = " << nEvents-nEventsRecovered<< endl;
	500	filcout << prefixDA << " : Nb of events used = " << nEvents << endl;
	501
	502	// Computing gain
	503	if(nIndex==nEntries)
dc928fbf	504	{
	505	muonGain->SetAliInit(nInit); // fnInit
	506	muonGain->SetAliEntries(nEntries); // fnEntries
	507	muonGain->SetAliNbpf1(nbpf1); // fnbpf1
1dd09aa3	508	muonGain->MakeGainStore(shuttleFile);
dc928fbf	509	}
dc928fbf	510
1dd09aa3	511	// Copying files to local DB folder defined by DAQ_DETDB_LOCAL
1dd09aa3	512	Char_t *dir;
dc928fbf	513	dir= getenv("DAQ_DETDB_LOCAL");
	514	unsigned int nLastVersions = 2;
	515	cout << "\n * Output files stored locally in " << dir << " (nb of previous versions = " << nLastVersions << ") *" << endl;
1dd09aa3	516	filcout << "\n * Output files stored locally in " << dir << " (nb of previous versions = " << nLastVersions << ") *" << endl;
	517
	518	filcout << endl;
	519	filcout << prefixDA << " : Root data file : " << muonGain->GetRootDataFileName() << endl;
	520	filcout << prefixDA << " : Output logfile : " << logOutputFile << endl;
	521	filcout << prefixDA << " : Gain Histo file : " << muonGain->GetHistoFileName() << endl;
	522	filcout << prefixDA << " : Gain file (to SHUTTLE) : " << shuttleFile << endl;
dc928fbf	523
	524	status = daqDA_localDB_storeFile(logOutputFile.Data(),nLastVersions);
	525	if(nIndex==nEntries)status = daqDA_localDB_storeFile(muonGain->GetRootDataFileName(),nLastVersions);
	526	status = daqDA_localDB_storeFile(muonGain->GetHistoFileName(),nLastVersions);
	527	status = daqDA_localDB_storeFile(shuttleFile.Data(),nLastVersions);
	528
	529	}
	530
dc928fbf	531	// ouput files
	532	cout << endl;
	533	cout << prefixDA << " : Root data file : " << muonGain->GetRootDataFileName() << endl;
	534	cout << prefixDA << " : Output logfile : " << logOutputFile << endl;
	535	cout << prefixDA << " : Gain Histo file : " << muonGain->GetHistoFileName() << endl;
	536	cout << prefixDA << " : Gain file (to SHUTTLE) : " << shuttleFile << endl;
1dd09aa3	537
dc928fbf	538
	539	// Transferring to OCDB via the SHUTTLE
	540	printf("\n *** STORE FILE in FES **** \n");
	541
	542	// be sure that env variable DAQDALIB_PATH is set in script file
	543	// gSystem->Setenv("DAQDALIB_PATH", "$DATE_SITE/infoLogger");
	544
	545	status = daqDA_FES_storeFile(shuttleFile.Data(),"GAINS");
	546	if (status)
	547	{
	548	printf(" Failed to export file : %d\n",status);
	549	}
	550	else printf(" %s successfully exported to FES \n",shuttleFile.Data());
	551
1dd09aa3	552
dc928fbf	553	filcout.close();
	554	timers.Stop();
	555	printf("\nExecution time : R:%7.2fs C:%7.2fs\n", timers.RealTime(), timers.CpuTime());
	556	return status;
1dd09aa3	557	}