[u/mrichter/AliRoot.git] / MUON / MUONcalib / 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),
fHistos(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(),
fHistoFile(0),
fTree(0)
{
/// 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),
fHistos(0),
fErrorBuspatchTable(0),
fManuBuspatchTable(0),
fManuBPoutofconfigTable(0),
fDate(0),
fFilcout(0),
fHistoFileName(),
fPedestalStore(0),
fIndex(-1),
fPrefixDA(),
fPrefixLDC(),
fHistoFile(0),
fTree(0)
{
/// 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 tree_charge=0;
  static Int_t warn=0;
  Int_t DDL= busPatchId/100+2560;
  /// Compute pedestals values
  AliMUONVCalibParam* ped = 
    static_cast<AliMUONVCalibParam*>(fPedestalStore ->FindObject(busPatchId, manuId));


  if(!tree_charge && fHistos==2)
    {
      fTree = new TTree("tc","Charge tree");
      fTree->Branch("bp",&busPatchId,"bp/I");
      fTree->Branch("manu",&manuId,",manu/I");
      fTree->Branch("channel",&channelId,",channel/I");
      fTree->Branch("DDL",&DDL,",DDL/I");
      fTree->Branch("charge",&charge,"charge/I");
      //      fTree->Branch("Pedestal",&Pedestal,"Pedestal/D");
      //      fTree->Branch("chargetrue",&chargeminusPed,"chargetrue/D");
      //      fTree->Branch("evt",&evt,"evt/I");
      tree_charge=1;
    }


  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.);

  if(fHistos==2) fTree->Fill();  

  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::CreateControlHistos()
{
// Create histo 
  fHistoFileName=Form("%s.root",fPrefixDA.Data());
  fHistoFile = new TFile(fHistoFileName,"RECREATE","MUON Tracking pedestals");
}
//______________________________________________________________________________
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;
  Int_t DDL;

// 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(fHistoFile);

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

  tree = new TTree("t","Pedestal tree");
  tree->Branch("DDL",&DDL,",DDL/I");
  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();
    DDL= busPatchId/100+2560;
    
    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();  
    }
  }
    
  fHistoFile->Write();  
  fHistoFile->Close(); 

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