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

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

#include "AliMUONPedestal.h"
#include "AliMUONErrorCounter.h"
#include "AliMUONVStore.h"
#include "AliMUON2DMap.h"
#include "AliMUONCalibParamND.h"
#include "AliMpConstants.h"
#include <TString.h>
#include <TTimeStamp.h>
#include <TMath.h>
#include <TTree.h>
#include <TFile.h>
#include <TH1F.h>
#include <THashList.h>
#include <Riostream.h>

#include <sstream>

//-----------------------------------------------------------------------------
/// \class AliMUONPedestal
///
/// Implementation of the pedestal computing
///
/// add
/// 
///
/// \author Alberto Baldisseri, JL Charvet (05/05/2009)
//-----------------------------------------------------------------------------

using std::ostream;
using std::ifstream;
using std::endl;
using std::cout;
using std::ios;
/// \cond CLASSIMP
ClassImp(AliMUONPedestal)
/// \endcond

//______________________________________________________________________________
AliMUONPedestal::AliMUONPedestal()
: TObject(),
//fN(0),
fNCurrentEvents(0),
fNEvthreshold(0),
fSorting(0),
fNEvents(0),
fRunNumber(0),
fNChannel(0),
fNManu(0),
fNManuConfig(0),
fConfig(1),
fStatusDA(0),
fErrorBuspatchTable(new AliMUON2DMap(kFALSE)),
fManuBuspatchTable(new AliMUON2DMap(kFALSE)),
fManuBPoutofconfigTable(new AliMUON2DMap(kFALSE)),
fDate(new TTimeStamp()),
fFilcout(0),
fHistoFileName(),
fPedestalStore(new AliMUON2DMap(kTRUE)),
fIndex(-1),
fPrefixDA(),
fPrefixLDC()
{
/// Default constructor
}
//______________________________________________________________________________
AliMUONPedestal::AliMUONPedestal(TRootIOCtor* /*dummy*/)
: TObject(),
//fN(0),
fNCurrentEvents(0),
fNEvthreshold(0),
fSorting(0),
fNEvents(0),
fRunNumber(0),
fNChannel(0),
fNManu(0),
fNManuConfig(0),
fConfig(1),
fStatusDA(0),
fErrorBuspatchTable(0),
fManuBuspatchTable(0),
fManuBPoutofconfigTable(0),
fDate(0),
fFilcout(0),
fHistoFileName(),
fPedestalStore(0),
fIndex(-1),
fPrefixDA(),
fPrefixLDC()
{
/// Root IO constructor
}

//______________________________________________________________________________
AliMUONPedestal::~AliMUONPedestal()
{
/// Destructor
  delete fErrorBuspatchTable;
  delete fManuBuspatchTable;
  delete fPedestalStore;
  delete fManuBPoutofconfigTable;
}

//______________________________________________________________________________
const char* 
AliMUONPedestal::GetHistoFileName() const
{
  /// Return the name of file we use to store histograms
  return fHistoFileName.Data();
}

//______________________________________________________________________________
void AliMUONPedestal::LoadConfig(const char* dbfile)
{
  /// Load MUONTRK configuration from ascii file "dbfile" (in DetDB)

  Int_t manuId;
  Int_t busPatchId;

  ifstream filein(dbfile,ios::in);
  
  while (!filein.eof())
    { 
      filein >> busPatchId >> manuId;

      AliMUONVCalibParam* ped = 
	static_cast<AliMUONVCalibParam*>(fPedestalStore ->FindObject(busPatchId, manuId));

      if (!ped) {
	fNManuConfig++;
	fNChannel+=64;
  ped = new AliMUONCalibParamND(2, AliMpConstants::ManuNofChannels(),busPatchId, manuId, -1.); // put default wise -1, not connected channel
	fPedestalStore ->Add(ped);  

        if ( ! fManuBuspatchTable->FindObject(busPatchId,manuId) )
	  {
	    // New (buspatch,manu)
	    AliMUONErrorCounter* manuCounter = new AliMUONErrorCounter(busPatchId,manuId);
	    fManuBuspatchTable->Add(manuCounter);
	  }
      }
    }
} 
//______________________________________________________________________________
void AliMUONPedestal::MakePed(Int_t busPatchId, Int_t manuId, Int_t channelId, Int_t charge)
{
  static Int_t warn=0;
  /// Compute pedestals values
  AliMUONVCalibParam* ped = 
    static_cast<AliMUONVCalibParam*>(fPedestalStore ->FindObject(busPatchId, manuId));

  if (!ped)   
    {
      if(fConfig) 
	{  // Fill out_of_config (buspatch,manu) table
	  if (!(static_cast<AliMUONErrorCounter*>(fManuBPoutofconfigTable->FindObject(busPatchId,manuId))))
	    fManuBPoutofconfigTable->Add(new AliMUONErrorCounter(busPatchId,manuId));
	  if(warn<10) cout << fPrefixLDC.Data() << " : !!! WARNING  : busPatchId = " << busPatchId << " manuId = " << manuId << " not in the Detector configuration " << endl;
	  else if(warn==10) cout << fPrefixLDC.Data() << " : !!! see .log file for an exhaustive list of (busPatchId, manuId) out of Detector configuration \n" << endl; 
	   warn++;
	   (*fFilcout) <<" !!! WARNING  : busPatchId = " << busPatchId << " manuId = " << manuId << " not in the Detector configuration " << endl; 
	}
      else {fNManu++;}
      fNChannel+=64;
      // put default wise -1, not connected channel
      ped = new AliMUONCalibParamND(2, AliMpConstants::ManuNofChannels(),busPatchId, manuId, -1.); 
      fPedestalStore ->Add(ped);  
    }

  // Initialization for the first value
  if (ped->ValueAsDouble(channelId, 0) == -1)  
    { 
      if(fConfig && channelId == 0){fNManu++;}
      ped->SetValueAsDouble(channelId, 0, 0.);
    }
  if (ped->ValueAsDouble(channelId, 1) == -1) ped->SetValueAsDouble(channelId, 1, 0.);

  Double_t pedMean  = ped->ValueAsDouble(channelId, 0) + (Double_t) charge;
  Double_t pedSigma = ped->ValueAsDouble(channelId, 1) + (Double_t) charge*charge;

  ped->SetValueAsDouble(channelId, 0, pedMean);
  ped->SetValueAsDouble(channelId, 1, pedSigma);

  AliMUONErrorCounter* manuCounter;
  if (!(manuCounter = static_cast<AliMUONErrorCounter*>(fManuBuspatchTable->FindObject(busPatchId,manuId))))
    {
      // New (buspatch,manu)
      manuCounter = new AliMUONErrorCounter(busPatchId,manuId);
      fManuBuspatchTable->Add(manuCounter);
    }
  else
    {
      // Existing buspatch
      manuCounter->Increment();
    }	
}
//______________________________________________________________________________
void AliMUONPedestal::Finalize()
{
  /// final polishing of the store
  
  Double_t pedMean;
  Double_t pedSigma;
  Double_t pedSigmalimit=0.5;
  Int_t busPatchId;
  Int_t manuId;
  Int_t channelId;
  Int_t status=0;

  // print in logfile
  if (fErrorBuspatchTable->GetSize())
    {
      cout<< "\n" << fPrefixLDC.Data() << " : See list of Buspatches with lower statistics (due to parity errors) in .log or .ped file "<<endl;
      (*fFilcout)<<"\nWarning: Buspatches with less statistics (due to parity errors)"<<endl;
      TIter nextParityError(fErrorBuspatchTable->CreateIterator());
      AliMUONErrorCounter* parityerror;
      while((parityerror = static_cast<AliMUONErrorCounter*>(nextParityError())))
	{
	  //	  cout<<"  bp "<<parityerror->BusPatch()<<": used events = "<<fNEvents-parityerror->Events()<<endl;
	  (*fFilcout)<<"  bp "<<parityerror->BusPatch()<<": used events = "<<fNEvents-parityerror->Events()<<endl;
	}
    }
  
  Int_t nADC4090=0;
  Int_t nADCmax=0;
  // iterator over pedestal
  TIter next(fPedestalStore ->CreateIterator());
  AliMUONVCalibParam* ped;

  while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
    {
      busPatchId              = ped->ID0();
      manuId                  = ped->ID1();
      if(manuId==0)
	{
	  cout << fPrefixLDC.Data() << " : Warning: ManuId = " << manuId << " !!! in  BP = " << busPatchId << endl;
	  (*fFilcout) << "Warning: ManuId = " << manuId << " !!! in  BP = " << busPatchId << endl;
	}
      Int_t eventCounter;
      // Correct the number of events for buspatch with errors
      AliMUONErrorCounter* errorCounter;
      if ((errorCounter = (AliMUONErrorCounter*)fErrorBuspatchTable->FindObject(busPatchId)))
	{
	  eventCounter = fNEvents - errorCounter->Events();
	}
      else
	{
	  eventCounter = fNEvents;
	}
      Int_t occupancy=0; // channel missing in raw data or read but rejected (case of parity error)
      // value of (buspatch, manu) occupancy
      AliMUONErrorCounter* manuCounter;
      manuCounter = static_cast<AliMUONErrorCounter*>(fManuBuspatchTable->FindObject(busPatchId,manuId));
      if(eventCounter>0)occupancy = manuCounter->Events()/64/eventCounter;
      if(occupancy>1)
	{
	  cout << fPrefixLDC.Data() << " : Warning: ManuId = " << manuId << " !!! in  BP = " << busPatchId << " occupancy (>1) = " << occupancy << endl;
	  (*fFilcout) << "Warning: ManuId = " << manuId << " !!! in  BP = " << busPatchId << " occupancy (>1) = " << occupancy <<endl;
	}

      for (channelId = 0; channelId < ped->Size() ; ++channelId) 
	{
	  pedMean  = ped->ValueAsDouble(channelId, 0);

	  if (pedMean >= 0) // connected channels
	    {
	      ped->SetValueAsDouble(channelId, 0, pedMean/(Double_t)eventCounter);
	      pedMean  = ped->ValueAsDouble(channelId, 0);
	      pedSigma = ped->ValueAsDouble(channelId, 1);
	      ped->SetValueAsDouble(channelId, 1, TMath::Sqrt(TMath::Abs(pedSigma/(Double_t)eventCounter - pedMean*pedMean)));

	      if(eventCounter < fNEvthreshold )
		{ nADCmax++; ped->SetValueAsDouble(channelId, 0, ADCMax());
		  ped->SetValueAsDouble(channelId, 1, ADCMax());}
	      if( ped->ValueAsDouble(channelId, 1) < pedSigmalimit )
		{ nADC4090++; ped->SetValueAsDouble(channelId, 0, ADCMax()-5);
		  ped->SetValueAsDouble(channelId, 1, ADCMax()-5);}
	      if(manuId == 0 || occupancy>1)
		{ nADCmax++; ped->SetValueAsDouble(channelId, 0, ADCMax());
		  ped->SetValueAsDouble(channelId, 1, ADCMax());
		  if(occupancy>1 && channelId==0)ped->SetValueAsDouble(channelId, 0, ADCMax()+occupancy);}
	    }
	  else
	    { nADCmax++; ped->SetValueAsDouble(channelId, 0, ADCMax());
	      ped->SetValueAsDouble(channelId, 1, ADCMax());}
	}
    }

  float frac1=0. , frac2=0. ;
  float frac_badped = 0.25 ; // maximal acceptable ratio of bad pedestals
  char* detail;
 
  if(fNChannel)
    {
      if(nADCmax>0)
	{ frac1 = float(nADCmax)/float(fNChannel); 
	  detail=Form("%s : Warning: Nb of Channels with bad Pedestal (Ped=4095) = %d over %d (%6.4f)",fPrefixLDC.Data(),nADCmax,fNChannel,frac1);
	  printf("%s\n",detail);
	  (*fFilcout) <<  detail << endl;}

      if(nADC4090>0)
	{ frac2 = 1.*nADC4090/fNChannel; 
	  detail=Form("%s : Warning: Nb of Channels with PedSigma<0.5 (Ped=4090) = %d over %d (%6.4f)",fPrefixLDC.Data(),nADC4090,fNChannel,frac2);
	  printf("%s\n",detail);
	  (*fFilcout) <<  detail << endl; }

      if (frac1+frac2 > frac_badped) { status=-1 ;
	detail=Form("\n%s !!! ERROR : fraction of Channels with Pedestal>=4090 = %6.4f (> %5.3f)  (status= %d) \n",fPrefixLDC.Data(),frac1+frac2,frac_badped,status);
	(*fFilcout) <<  detail << endl; printf("%s",detail) ;}
    }
  else { status= -1;
    detail=Form("\n%s !!! ERROR : Nb good channel = 0 (all pedestals are forced to 4095) !!! (status= %d)\n",fPrefixLDC.Data(),status); 
    cout << detail; (*fFilcout) << detail ;}
   
  SetStatusDA(status);
}
//______________________________________________________________________________
void AliMUONPedestal::MakeASCIIoutput(ostream& out) const
{
/// put pedestal store in the output stream

  out<<"//===========================================================================" << endl;
  out<<"//                 Pedestal file calculated by "<< fPrefixDA.Data() << endl;
  out<<"//===========================================================================" << endl;
  out<<"//     * Run           : " << fRunNumber << endl; 
  out<<"//     * Date          : " << fDate->AsString("l") <<endl;
  out<<"//     * Statictics    : " << fNEvents << endl;
  if(fConfig)
    out<<"//     * Nb of MANUS   : " << fNManuConfig << " read in the Det. config. " << endl;
  out<<"//     * Nb of MANUS   : " << fNManu << " read in raw data " << endl;
  out<<"//     * Nb of MANUS   : " << fNChannel/64 << " written in pedestal file " << endl;
  out<<"//     * Nb of channels: " << fNChannel << endl;
  out<<"//"<<endl;
  out<<"//     * Below " << fNEvthreshold << " events=> Ped.&sig.=4095" << endl;
  out<<"//     * PedSigma < 0.5 => Ped.&sig.=4090" << endl;

  if (fErrorBuspatchTable->GetSize())
    {
      out<<"//"<<endl;
      out<<"//    * Buspatches with less statistics (due to parity errors)"<<endl;
      TIter next(fErrorBuspatchTable->CreateIterator());
      AliMUONErrorCounter* parityerror;
      while((parityerror = static_cast<AliMUONErrorCounter*>(next())))
	{
	  if(fNEvents-parityerror->Events()>fNEvthreshold)
	    { out<<"//      BusPatch = "<<parityerror->BusPatch()<<"\t Nevents used = "<<fNEvents-parityerror->Events()<<endl; }
	  else
	    { out<<"//      BusPatch = "<<parityerror->BusPatch()<<"\t Nevents used = "<<fNEvents-parityerror->Events()<< " (Ped.&sig.=4095)" << endl; }
	}
    }  
  Int_t writitle=0;
  Int_t occupancy=1;
  if(occupancy)
    {
      TIter next(fPedestalStore ->CreateIterator());
      AliMUONVCalibParam* ped;
      while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
	{
	  Int_t busPatchId = ped->ID0();
	  Int_t manuId = ped->ID1();
	  Double_t pedMean  = ped->ValueAsDouble(0, 0); // check pedestal value for channelId=0

	  if(pedMean>ADCMax()) 
	    {
	      writitle++;
	      if(writitle==1){ 
		out<<"//"<<endl;
		out<<"//    * Puzzling (Buspatch,Manu) read in raw data ? (Ped.&sig.=4095)"<<endl;}
	      occupancy=TMath::Nint(pedMean-ADCMax());
	      ped->SetValueAsDouble(0, 0, ADCMax());
	      out<<"//      BusPatch = "<< busPatchId <<"\t ManuId =  "<< manuId << "\t occupancy = " << occupancy  << endl;
	    }

	  if (manuId==0 || (fConfig && static_cast<AliMUONErrorCounter*>(fManuBPoutofconfigTable->FindObject(busPatchId,manuId))))
	    {
	      writitle++;
	      if(writitle==1){ 
		out<<"//"<<endl;
		out<<"//    * Puzzling (Buspatch,Manu) read in raw data ? (Ped.&sig.=4095)"<<endl;}
	      out<<"//      BusPatch = "<< busPatchId <<"\t ManuId =  "<< manuId << "\t missing in the mapping" << endl;
	    }
	}
    }
  out<<"//"<<endl;
  out<<"//---------------------------------------------------------------------------" << endl;
  out<<"//---------------------------------------------------------------------------" << endl;
  out<<"//      BP     MANU     CH.      MEAN    SIGMA"<<endl;
  out<<"//---------------------------------------------------------------------------" << endl;

  TIter next(fPedestalStore->CreateIterator());
  AliMUONVCalibParam* ped;  

  // Sorting 
  if  (fSorting)
    {
      printf("%s : ..... sorting pedestal values .....\n",fPrefixLDC.Data());
      THashList pedtable(100,2);
      while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
      {
        pedtable.Add(ped);
      }
      pedtable.Sort();
      //      iterator over sorted pedestal
      TIter nextSorted(&pedtable);
      while ( (ped = (AliMUONVCalibParam*)(nextSorted()) ) )
      {
        Int_t busPatchId = ped->ID0();
        Int_t manuId = ped->ID1();
        for ( Int_t channelId = 0; channelId < ped->Size(); ++channelId ) 
        {
          Double_t pedMean  = ped->ValueAsDouble(channelId, 0);
          Double_t pedSigma = ped->ValueAsDouble(channelId, 1);
          out << "\t" << busPatchId << "\t" << manuId <<"\t"<< channelId << "\t" << pedMean <<"\t"<< pedSigma << endl;
        }
      }
    }
  else
    {
      while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
	{
	  Int_t busPatchId = ped->ID0();
	  Int_t manuId = ped->ID1();
	  for ( Int_t channelId = 0; channelId < ped->Size(); ++channelId ) 
	    {
	      Double_t pedMean  = ped->ValueAsDouble(channelId, 0);
	      Double_t pedSigma = ped->ValueAsDouble(channelId, 1);
	      out << "\t" << busPatchId << "\t" << manuId <<"\t"<< channelId << "\t" << pedMean <<"\t"<< pedSigma << endl;
	    }
	}
    }
}

//______________________________________________________________________________
void AliMUONPedestal::MakeControlHistos()
{
  /// Create control histograms
  if (fIndex>=0) return; // Pedestal run (fIndex=-1)

  Double_t pedMean;
  Double_t pedSigma;
  Double_t evt;
  Int_t busPatchId;
  Int_t manuId;
  Int_t channelId;

// histo
  TFile*  histoFile = 0;
  TTree* tree = 0;
  TH1F* pedMeanHisto = 0;
  TH1F* pedSigmaHisto = 0;
    
  fHistoFileName=Form("%s.root",fPrefixDA.Data());
  histoFile = new TFile(fHistoFileName,"RECREATE","MUON Tracking pedestals");

  Int_t nx = ADCMax()+1;
  Int_t xmin = 0;
  Int_t xmax = ADCMax(); 
  pedMeanHisto = new TH1F("pedmean_allch","Pedestal mean all channels",nx,xmin,xmax);
  pedMeanHisto->SetDirectory(histoFile);

  nx = 201;
  xmin = 0;
  xmax = 200; 
  pedSigmaHisto = new TH1F("pedsigma_allch","Pedestal sigma all channels",nx,xmin,xmax);
  pedSigmaHisto->SetDirectory(histoFile);

  tree = new TTree("t","Pedestal tree");
  tree->Branch("bp",&busPatchId,"bp/I");
  tree->Branch("manu",&manuId,",manu/I");
  tree->Branch("channel",&channelId,",channel/I");
  tree->Branch("pedMean",&pedMean,",pedMean/D");
  tree->Branch("pedSigma",&pedSigma,",pedSigma/D");
  tree->Branch("nevt",&evt,",evt/D");

  // iterator over pedestal
  TIter next(fPedestalStore ->CreateIterator());
  AliMUONVCalibParam* ped;
  AliMUONErrorCounter* manuCounter;
  
  while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
  {
    busPatchId = ped->ID0();
    manuId = ped->ID1();
    
    for ( channelId = 0; channelId < ped->Size(); ++channelId ) 
    {
      pedMean  = ped->ValueAsDouble(channelId, 0);
      pedSigma = ped->ValueAsDouble(channelId, 1);
      manuCounter = static_cast<AliMUONErrorCounter*>(fManuBuspatchTable->FindObject(busPatchId,manuId));
      evt = manuCounter->Events()/64;
          
      pedMeanHisto->Fill(pedMean);
      pedSigmaHisto->Fill(pedSigma);
      tree->Fill();  
    }
  }
    
  histoFile->Write();  
  histoFile->Close(); 

}
Commit	Line	Data
5253c153	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	#include "AliMUONPedestal.h"
	19	#include "AliMUONErrorCounter.h"
	20	#include "AliMUONVStore.h"
	21	#include "AliMUON2DMap.h"
	22	#include "AliMUONCalibParamND.h"
f3cfa63e	23	#include "AliMpConstants.h"
5253c153	24	#include <TString.h>
5253c153	25	#include <TTimeStamp.h>
	26	#include <TMath.h>
	27	#include <TTree.h>
	28	#include <TFile.h>
	29	#include <TH1F.h>
ff496840	30	#include <THashList.h>
5253c153	31	#include <Riostream.h>
	32
	33	#include <sstream>
	34
	35	//-----------------------------------------------------------------------------
	36	/// \class AliMUONPedestal
	37	///
	38	/// Implementation of the pedestal computing
	39	///
	40	/// add
	41	///
	42	///
	43	/// \author Alberto Baldisseri, JL Charvet (05/05/2009)
	44	//-----------------------------------------------------------------------------
	45
b80faac0	46	using std::ostream;
	47	using std::ifstream;
	48	using std::endl;
	49	using std::cout;
	50	using std::ios;
5253c153	51	/// \cond CLASSIMP
	52	ClassImp(AliMUONPedestal)
	53	/// \endcond
	54
	55	//______________________________________________________________________________
	56	AliMUONPedestal::AliMUONPedestal()
	57	: TObject(),
f3cfa63e	58	//fN(0),
f9e83f1d	59	fNCurrentEvents(0),
f9e83f1d	60	fNEvthreshold(0),
ff496840	61	fSorting(0),
5253c153	62	fNEvents(0),
	63	fRunNumber(0),
	64	fNChannel(0),
	65	fNManu(0),
f3cfa63e	66	fNManuConfig(0),
d9e6dc5e	67	fConfig(1),
9c2d4e89	68	fStatusDA(0),
a54854fc	69	fErrorBuspatchTable(new AliMUON2DMap(kFALSE)),
a54854fc	70	fManuBuspatchTable(new AliMUON2DMap(kFALSE)),
1b15b395	71	fManuBPoutofconfigTable(new AliMUON2DMap(kFALSE)),
5253c153	72	fDate(new TTimeStamp()),
5253c153	73	fFilcout(0),
f3cfa63e	74	fHistoFileName(),
a54854fc	75	fPedestalStore(new AliMUON2DMap(kTRUE)),
f3cfa63e	76	fIndex(-1),
9c2d4e89	77	fPrefixDA(),
9c2d4e89	78	fPrefixLDC()
5253c153	79	{
5253c153	80	/// Default constructor
5253c153	81	}
71432c59	82	//______________________________________________________________________________
	83	AliMUONPedestal::AliMUONPedestal(TRootIOCtor* /dummy/)
	84	: TObject(),
	85	//fN(0),
	86	fNCurrentEvents(0),
	87	fNEvthreshold(0),
	88	fSorting(0),
	89	fNEvents(0),
	90	fRunNumber(0),
	91	fNChannel(0),
	92	fNManu(0),
	93	fNManuConfig(0),
	94	fConfig(1),
9c2d4e89	95	fStatusDA(0),
71432c59	96	fErrorBuspatchTable(0),
	97	fManuBuspatchTable(0),
	98	fManuBPoutofconfigTable(0),
	99	fDate(0),
	100	fFilcout(0),
	101	fHistoFileName(),
	102	fPedestalStore(0),
	103	fIndex(-1),
9c2d4e89	104	fPrefixDA(),
9c2d4e89	105	fPrefixLDC()
71432c59	106	{
	107	/// Root IO constructor
	108	}
	109
5253c153	110	//______________________________________________________________________________
	111	AliMUONPedestal::~AliMUONPedestal()
	112	{
	113	/// Destructor
a54854fc	114	delete fErrorBuspatchTable;
a54854fc	115	delete fManuBuspatchTable;
dbbb2c64	116	delete fPedestalStore;
1b15b395	117	delete fManuBPoutofconfigTable;
5253c153	118	}
5253c153	119
1b15b395	120	//______________________________________________________________________________
f3cfa63e	121	const char*
	122	AliMUONPedestal::GetHistoFileName() const
	123	{
	124	/// Return the name of file we use to store histograms
	125	return fHistoFileName.Data();
	126	}
	127
	128	//______________________________________________________________________________
	129	void AliMUONPedestal::LoadConfig(const char* dbfile)
1b15b395	130	{
	131	/// Load MUONTRK configuration from ascii file "dbfile" (in DetDB)
	132
	133	Int_t manuId;
	134	Int_t busPatchId;
	135
f3cfa63e	136	ifstream filein(dbfile,ios::in);
f3cfa63e	137
f3cfa63e	138	while (!filein.eof())
1b15b395	139	{
f3cfa63e	140	filein >> busPatchId >> manuId;
1b15b395	141
1b15b395	142	AliMUONVCalibParam* ped =
	143	static_cast<AliMUONVCalibParam*>(fPedestalStore ->FindObject(busPatchId, manuId));
	144
	145	if (!ped) {
f3cfa63e	146	fNManuConfig++;
1b15b395	147	fNChannel+=64;
f3cfa63e	148	ped = new AliMUONCalibParamND(2, AliMpConstants::ManuNofChannels(),busPatchId, manuId, -1.); // put default wise -1, not connected channel
1b15b395	149	fPedestalStore ->Add(ped);
1b15b395	150
065a2d7d	151	if ( ! fManuBuspatchTable->FindObject(busPatchId,manuId) )
1b15b395	152	{
1b15b395	153	// New (buspatch,manu)
065a2d7d	154	AliMUONErrorCounter* manuCounter = new AliMUONErrorCounter(busPatchId,manuId);
1b15b395	155	fManuBuspatchTable->Add(manuCounter);
	156	}
	157	}
	158	}
	159	}
5253c153	160	//______________________________________________________________________________
	161	void AliMUONPedestal::MakePed(Int_t busPatchId, Int_t manuId, Int_t channelId, Int_t charge)
	162	{
b9c33bca	163	static Int_t warn=0;
5253c153	164	/// Compute pedestals values
	165	AliMUONVCalibParam* ped =
	166	static_cast<AliMUONVCalibParam*>(fPedestalStore ->FindObject(busPatchId, manuId));
	167
1b15b395	168	if (!ped)
	169	{
	170	if(fConfig)
	171	{ // Fill out_of_config (buspatch,manu) table
	172	if (!(static_cast<AliMUONErrorCounter*>(fManuBPoutofconfigTable->FindObject(busPatchId,manuId))))
	173	fManuBPoutofconfigTable->Add(new AliMUONErrorCounter(busPatchId,manuId));
9c2d4e89	174	if(warn<10) cout << fPrefixLDC.Data() << " : !!! WARNING : busPatchId = " << busPatchId << " manuId = " << manuId << " not in the Detector configuration " << endl;
9c2d4e89	175	else if(warn==10) cout << fPrefixLDC.Data() << " : !!! see .log file for an exhaustive list of (busPatchId, manuId) out of Detector configuration \n" << endl;
b9c33bca	176	warn++;
9c2d4e89	177	(*fFilcout) <<" !!! WARNING : busPatchId = " << busPatchId << " manuId = " << manuId << " not in the Detector configuration " << endl;
1b15b395	178	}
	179	else {fNManu++;}
	180	fNChannel+=64;
	181	// put default wise -1, not connected channel
f3cfa63e	182	ped = new AliMUONCalibParamND(2, AliMpConstants::ManuNofChannels(),busPatchId, manuId, -1.);
1b15b395	183	fPedestalStore ->Add(ped);
1b15b395	184	}
5253c153	185
5253c153	186	// Initialization for the first value
1b15b395	187	if (ped->ValueAsDouble(channelId, 0) == -1)
	188	{
	189	if(fConfig && channelId == 0){fNManu++;}
	190	ped->SetValueAsDouble(channelId, 0, 0.);
	191	}
5253c153	192	if (ped->ValueAsDouble(channelId, 1) == -1) ped->SetValueAsDouble(channelId, 1, 0.);
	193
	194	Double_t pedMean = ped->ValueAsDouble(channelId, 0) + (Double_t) charge;
	195	Double_t pedSigma = ped->ValueAsDouble(channelId, 1) + (Double_t) charge*charge;
	196
	197	ped->SetValueAsDouble(channelId, 0, pedMean);
	198	ped->SetValueAsDouble(channelId, 1, pedSigma);
5431405e	199
5431405e	200	AliMUONErrorCounter* manuCounter;
a54854fc	201	if (!(manuCounter = static_cast<AliMUONErrorCounter*>(fManuBuspatchTable->FindObject(busPatchId,manuId))))
5431405e	202	{
	203	// New (buspatch,manu)
	204	manuCounter = new AliMUONErrorCounter(busPatchId,manuId);
a54854fc	205	fManuBuspatchTable->Add(manuCounter);
5431405e	206	}
	207	else
	208	{
	209	// Existing buspatch
	210	manuCounter->Increment();
	211	}
5253c153	212	}
5253c153	213	//______________________________________________________________________________
dbbb2c64	214	void AliMUONPedestal::Finalize()
5253c153	215	{
f3cfa63e	216	/// final polishing of the store
f3cfa63e	217
5253c153	218	Double_t pedMean;
5253c153	219	Double_t pedSigma;
f9e83f1d	220	Double_t pedSigmalimit=0.5;
5253c153	221	Int_t busPatchId;
	222	Int_t manuId;
	223	Int_t channelId;
9c2d4e89	224	Int_t status=0;
5253c153	225
5253c153	226	// print in logfile
5253c153	227	if (fErrorBuspatchTable->GetSize())
5253c153	228	{
9c2d4e89	229	cout<< "\n" << fPrefixLDC.Data() << " : See list of Buspatches with lower statistics (due to parity errors) in .log or .ped file "<<endl;
f9e83f1d	230	(*fFilcout)<<"\nWarning: Buspatches with less statistics (due to parity errors)"<<endl;
1b15b395	231	TIter nextParityError(fErrorBuspatchTable->CreateIterator());
	232	AliMUONErrorCounter* parityerror;
	233	while((parityerror = static_cast<AliMUONErrorCounter*>(nextParityError())))
	234	{
9c2d4e89	235	// cout<<" bp "<<parityerror->BusPatch()<<": used events = "<<fNEvents-parityerror->Events()<<endl;
9c2d4e89	236	(*fFilcout)<<" bp "<<parityerror->BusPatch()<<": used events = "<<fNEvents-parityerror->Events()<<endl;
1b15b395	237	}
5253c153	238	}
9c2d4e89	239
f9e83f1d	240	Int_t nADC4090=0;
f9e83f1d	241	Int_t nADCmax=0;
1b15b395	242	// iterator over pedestal
5253c153	243	TIter next(fPedestalStore ->CreateIterator());
	244	AliMUONVCalibParam* ped;
	245
	246	while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
dbbb2c64	247	{
1b15b395	248	busPatchId = ped->ID0();
	249	manuId = ped->ID1();
	250	if(manuId==0)
	251	{
9c2d4e89	252	cout << fPrefixLDC.Data() << " : Warning: ManuId = " << manuId << " !!! in BP = " << busPatchId << endl;
f9e83f1d	253	(*fFilcout) << "Warning: ManuId = " << manuId << " !!! in BP = " << busPatchId << endl;
1b15b395	254	}
	255	Int_t eventCounter;
	256	// Correct the number of events for buspatch with errors
	257	AliMUONErrorCounter* errorCounter;
	258	if ((errorCounter = (AliMUONErrorCounter*)fErrorBuspatchTable->FindObject(busPatchId)))
	259	{
	260	eventCounter = fNEvents - errorCounter->Events();
	261	}
	262	else
	263	{
	264	eventCounter = fNEvents;
	265	}
1b15b395	266	Int_t occupancy=0; // channel missing in raw data or read but rejected (case of parity error)
	267	// value of (buspatch, manu) occupancy
	268	AliMUONErrorCounter* manuCounter;
	269	manuCounter = static_cast<AliMUONErrorCounter*>(fManuBuspatchTable->FindObject(busPatchId,manuId));
	270	if(eventCounter>0)occupancy = manuCounter->Events()/64/eventCounter;
	271	if(occupancy>1)
	272	{
9c2d4e89	273	cout << fPrefixLDC.Data() << " : Warning: ManuId = " << manuId << " !!! in BP = " << busPatchId << " occupancy (>1) = " << occupancy << endl;
f9e83f1d	274	(*fFilcout) << "Warning: ManuId = " << manuId << " !!! in BP = " << busPatchId << " occupancy (>1) = " << occupancy <<endl;
1b15b395	275	}
5253c153	276
1b15b395	277	for (channelId = 0; channelId < ped->Size() ; ++channelId)
1b15b395	278	{
9c2d4e89	279	pedMean = ped->ValueAsDouble(channelId, 0);
1b15b395	280
f9e83f1d	281	if (pedMean >= 0) // connected channels
1b15b395	282	{
	283	ped->SetValueAsDouble(channelId, 0, pedMean/(Double_t)eventCounter);
	284	pedMean = ped->ValueAsDouble(channelId, 0);
	285	pedSigma = ped->ValueAsDouble(channelId, 1);
9c2d4e89	286	ped->SetValueAsDouble(channelId, 1, TMath::Sqrt(TMath::Abs(pedSigma/(Double_t)eventCounter - pedMean*pedMean)));
f9e83f1d	287
	288	if(eventCounter < fNEvthreshold )
	289	{ nADCmax++; ped->SetValueAsDouble(channelId, 0, ADCMax());
	290	ped->SetValueAsDouble(channelId, 1, ADCMax());}
	291	if( ped->ValueAsDouble(channelId, 1) < pedSigmalimit )
	292	{ nADC4090++; ped->SetValueAsDouble(channelId, 0, ADCMax()-5);
	293	ped->SetValueAsDouble(channelId, 1, ADCMax()-5);}
	294	if(manuId == 0 \|\| occupancy>1)
	295	{ nADCmax++; ped->SetValueAsDouble(channelId, 0, ADCMax());
f3cfa63e	296	ped->SetValueAsDouble(channelId, 1, ADCMax());
f9e83f1d	297	if(occupancy>1 && channelId==0)ped->SetValueAsDouble(channelId, 0, ADCMax()+occupancy);}
1b15b395	298	}
1b15b395	299	else
f9e83f1d	300	{ nADCmax++; ped->SetValueAsDouble(channelId, 0, ADCMax());
f9e83f1d	301	ped->SetValueAsDouble(channelId, 1, ADCMax());}
1b15b395	302	}
dbbb2c64	303	}
9c2d4e89	304
	305	float frac1=0. , frac2=0. ;
	306	float frac_badped = 0.25 ; // maximal acceptable ratio of bad pedestals
	307	char* detail;
	308
	309	if(fNChannel)
	310	{
	311	if(nADCmax>0)
	312	{ frac1 = float(nADCmax)/float(fNChannel);
	313	detail=Form("%s : Warning: Nb of Channels with bad Pedestal (Ped=4095) = %d over %d (%6.4f)",fPrefixLDC.Data(),nADCmax,fNChannel,frac1);
	314	printf("%s\n",detail);
	315	(*fFilcout) << detail << endl;}
	316
	317	if(nADC4090>0)
	318	{ frac2 = 1.*nADC4090/fNChannel;
	319	detail=Form("%s : Warning: Nb of Channels with PedSigma<0.5 (Ped=4090) = %d over %d (%6.4f)",fPrefixLDC.Data(),nADC4090,fNChannel,frac2);
	320	printf("%s\n",detail);
	321	(*fFilcout) << detail << endl; }
	322
	323	if (frac1+frac2 > frac_badped) { status=-1 ;
	324	detail=Form("\n%s !!! ERROR : fraction of Channels with Pedestal>=4090 = %6.4f (> %5.3f) (status= %d) \n",fPrefixLDC.Data(),frac1+frac2,frac_badped,status);
	325	(*fFilcout) << detail << endl; printf("%s",detail) ;}
	326	}
	327	else { status= -1;
	328	detail=Form("\n%s !!! ERROR : Nb good channel = 0 (all pedestals are forced to 4095) !!! (status= %d)\n",fPrefixLDC.Data(),status);
	329	cout << detail; (*fFilcout) << detail ;}
	330
	331	SetStatusDA(status);
dbbb2c64	332	}
dbbb2c64	333	//______________________________________________________________________________
	334	void AliMUONPedestal::MakeASCIIoutput(ostream& out) const
	335	{
9c2d4e89	336	/// put pedestal store in the output stream
dbbb2c64	337
dbbb2c64	338	out<<"//===========================================================================" << endl;
f3cfa63e	339	out<<"// Pedestal file calculated by "<< fPrefixDA.Data() << endl;
dbbb2c64	340	out<<"//===========================================================================" << endl;
f9e83f1d	341	out<<"// * Run : " << fRunNumber << endl;
	342	out<<"// * Date : " << fDate->AsString("l") <<endl;
	343	out<<"// * Statictics : " << fNEvents << endl;
1b15b395	344	if(fConfig)
f9e83f1d	345	out<<"// * Nb of MANUS : " << fNManuConfig << " read in the Det. config. " << endl;
	346	out<<"// * Nb of MANUS : " << fNManu << " read in raw data " << endl;
	347	out<<"// * Nb of MANUS : " << fNChannel/64 << " written in pedestal file " << endl;
	348	out<<"// * Nb of channels: " << fNChannel << endl;
	349	out<<"//"<<endl;
	350	out<<"// * Below " << fNEvthreshold << " events=> Ped.&sig.=4095" << endl;
	351	out<<"// * PedSigma < 0.5 => Ped.&sig.=4090" << endl;
	352
dbbb2c64	353	if (fErrorBuspatchTable->GetSize())
dbbb2c64	354	{
1b15b395	355	out<<"//"<<endl;
	356	out<<"// * Buspatches with less statistics (due to parity errors)"<<endl;
	357	TIter next(fErrorBuspatchTable->CreateIterator());
	358	AliMUONErrorCounter* parityerror;
	359	while((parityerror = static_cast<AliMUONErrorCounter*>(next())))
	360	{
f9e83f1d	361	if(fNEvents-parityerror->Events()>fNEvthreshold)
	362	{ out<<"// BusPatch = "<<parityerror->BusPatch()<<"\t Nevents used = "<<fNEvents-parityerror->Events()<<endl; }
	363	else
	364	{ out<<"// BusPatch = "<<parityerror->BusPatch()<<"\t Nevents used = "<<fNEvents-parityerror->Events()<< " (Ped.&sig.=4095)" << endl; }
1b15b395	365	}
1b15b395	366	}
f3cfa63e	367	Int_t writitle=0;
	368	Int_t occupancy=1;
	369	if(occupancy)
1b15b395	370	{
f3cfa63e	371	TIter next(fPedestalStore ->CreateIterator());
	372	AliMUONVCalibParam* ped;
	373	while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
1b15b395	374	{
f3cfa63e	375	Int_t busPatchId = ped->ID0();
	376	Int_t manuId = ped->ID1();
	377	Double_t pedMean = ped->ValueAsDouble(0, 0); // check pedestal value for channelId=0
1b15b395	378
f3cfa63e	379	if(pedMean>ADCMax())
	380	{
	381	writitle++;
	382	if(writitle==1){
	383	out<<"//"<<endl;
f9e83f1d	384	out<<"// * Puzzling (Buspatch,Manu) read in raw data ? (Ped.&sig.=4095)"<<endl;}
1ccd531d	385	occupancy=TMath::Nint(pedMean-ADCMax());
f3cfa63e	386	ped->SetValueAsDouble(0, 0, ADCMax());
f9e83f1d	387	out<<"// BusPatch = "<< busPatchId <<"\t ManuId = "<< manuId << "\t occupancy = " << occupancy << endl;
f3cfa63e	388	}
	389
	390	if (manuId==0 \|\| (fConfig && static_cast<AliMUONErrorCounter*>(fManuBPoutofconfigTable->FindObject(busPatchId,manuId))))
	391	{
	392	writitle++;
	393	if(writitle==1){
	394	out<<"//"<<endl;
f9e83f1d	395	out<<"// * Puzzling (Buspatch,Manu) read in raw data ? (Ped.&sig.=4095)"<<endl;}
f3cfa63e	396	out<<"// BusPatch = "<< busPatchId <<"\t ManuId = "<< manuId << "\t missing in the mapping" << endl;
f3cfa63e	397	}
1b15b395	398	}
dbbb2c64	399	}
dbbb2c64	400	out<<"//"<<endl;
	401	out<<"//---------------------------------------------------------------------------" << endl;
	402	out<<"//---------------------------------------------------------------------------" << endl;
	403	out<<"// BP MANU CH. MEAN SIGMA"<<endl;
	404	out<<"//---------------------------------------------------------------------------" << endl;
5253c153	405
ff496840	406	TIter next(fPedestalStore->CreateIterator());
ff496840	407	AliMUONVCalibParam* ped;
1b15b395	408
ff496840	409	// Sorting
	410	if (fSorting)
	411	{
9c2d4e89	412	printf("%s : ..... sorting pedestal values .....\n",fPrefixLDC.Data());
ff496840	413	THashList pedtable(100,2);
	414	while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
	415	{
	416	pedtable.Add(ped);
	417	}
	418	pedtable.Sort();
	419	// iterator over sorted pedestal
	420	TIter nextSorted(&pedtable);
	421	while ( (ped = (AliMUONVCalibParam*)(nextSorted()) ) )
	422	{
	423	Int_t busPatchId = ped->ID0();
	424	Int_t manuId = ped->ID1();
	425	for ( Int_t channelId = 0; channelId < ped->Size(); ++channelId )
	426	{
	427	Double_t pedMean = ped->ValueAsDouble(channelId, 0);
	428	Double_t pedSigma = ped->ValueAsDouble(channelId, 1);
	429	out << "\t" << busPatchId << "\t" << manuId <<"\t"<< channelId << "\t" << pedMean <<"\t"<< pedSigma << endl;
	430	}
	431	}
	432	}
	433	else
	434	{
	435	while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
1b15b395	436	{
ff496840	437	Int_t busPatchId = ped->ID0();
	438	Int_t manuId = ped->ID1();
	439	for ( Int_t channelId = 0; channelId < ped->Size(); ++channelId )
	440	{
	441	Double_t pedMean = ped->ValueAsDouble(channelId, 0);
	442	Double_t pedSigma = ped->ValueAsDouble(channelId, 1);
	443	out << "\t" << busPatchId << "\t" << manuId <<"\t"<< channelId << "\t" << pedMean <<"\t"<< pedSigma << endl;
	444	}
1b15b395	445	}
5253c153	446	}
dbbb2c64	447	}
	448
	449	//______________________________________________________________________________
	450	void AliMUONPedestal::MakeControlHistos()
	451	{
f3cfa63e	452	/// Create control histograms
dbbb2c64	453	if (fIndex>=0) return; // Pedestal run (fIndex=-1)
	454
	455	Double_t pedMean;
	456	Double_t pedSigma;
f9e83f1d	457	Double_t evt;
dbbb2c64	458	Int_t busPatchId;
	459	Int_t manuId;
	460	Int_t channelId;
	461
	462	// histo
	463	TFile* histoFile = 0;
	464	TTree* tree = 0;
	465	TH1F* pedMeanHisto = 0;
	466	TH1F* pedSigmaHisto = 0;
	467
f3cfa63e	468	fHistoFileName=Form("%s.root",fPrefixDA.Data());
dbbb2c64	469	histoFile = new TFile(fHistoFileName,"RECREATE","MUON Tracking pedestals");
dbbb2c64	470
f3cfa63e	471	Int_t nx = ADCMax()+1;
dbbb2c64	472	Int_t xmin = 0;
f3cfa63e	473	Int_t xmax = ADCMax();
eae882b9	474	pedMeanHisto = new TH1F("pedmean_allch","Pedestal mean all channels",nx,xmin,xmax);
dbbb2c64	475	pedMeanHisto->SetDirectory(histoFile);
dbbb2c64	476
dbbb2c64	477	nx = 201;
	478	xmin = 0;
	479	xmax = 200;
eae882b9	480	pedSigmaHisto = new TH1F("pedsigma_allch","Pedestal sigma all channels",nx,xmin,xmax);
dbbb2c64	481	pedSigmaHisto->SetDirectory(histoFile);
	482
	483	tree = new TTree("t","Pedestal tree");
	484	tree->Branch("bp",&busPatchId,"bp/I");
	485	tree->Branch("manu",&manuId,",manu/I");
	486	tree->Branch("channel",&channelId,",channel/I");
	487	tree->Branch("pedMean",&pedMean,",pedMean/D");
	488	tree->Branch("pedSigma",&pedSigma,",pedSigma/D");
f9e83f1d	489	tree->Branch("nevt",&evt,",evt/D");
dbbb2c64	490
	491	// iterator over pedestal
	492	TIter next(fPedestalStore ->CreateIterator());
	493	AliMUONVCalibParam* ped;
f9e83f1d	494	AliMUONErrorCounter* manuCounter;
dbbb2c64	495
	496	while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
	497	{
	498	busPatchId = ped->ID0();
	499	manuId = ped->ID1();
	500
	501	for ( channelId = 0; channelId < ped->Size(); ++channelId )
	502	{
	503	pedMean = ped->ValueAsDouble(channelId, 0);
	504	pedSigma = ped->ValueAsDouble(channelId, 1);
f9e83f1d	505	manuCounter = static_cast<AliMUONErrorCounter*>(fManuBuspatchTable->FindObject(busPatchId,manuId));
f9e83f1d	506	evt = manuCounter->Events()/64;
dbbb2c64	507
	508	pedMeanHisto->Fill(pedMean);
	509	pedSigmaHisto->Fill(pedSigma);
	510	tree->Fill();
	511	}
5253c153	512	}
dbbb2c64	513
	514	histoFile->Write();
	515	histoFile->Close();
	516
4af2c34a	517	}