[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)
//-----------------------------------------------------------------------------

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