[u/mrichter/AliRoot.git] / PMD / AliPMDCalibrator.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.		    *
// **************************************************************************/
//
//////////////////////////////////////////////////////////////////////////////

#include "TF1.h"
#include "TFile.h"
#include "TObjString.h"
#include "TROOT.h"
#include "TClonesArray.h"
#include "TH1F.h"
#include "TObjArray.h"
#include "TMath.h"

// --- Standard library ---

// --- AliRoot header files ---
#include "AliLog.h"
#include "AliRawReaderFile.h"
#include "AliPMDCalibrator.h"
#include "AliRawReaderDate.h"
#include "AliPMDRawStream.h"
#include "AliPMDCalibData.h"
#include "AliPMDddldata.h"
#include "AliCDBManager.h"
#include "AliCDBId.h"
#include "AliCDBMetaData.h"
#include "AliDAQ.h"

//#include "AliPMDCleanNoise.h"
//#include "AliPMDCleaner.h"

#include "AliPMDPedestal.h"
#include "AliCDBManager.h"
#include "AliCDBEntry.h"

ClassImp(AliPMDCalibrator)


AliPMDCalibrator::AliPMDCalibrator():
  fCalibGain(new AliPMDCalibData()), fCalibPed(new AliPMDPedestal())
{
  // Standard Constructor
  for(Int_t idet = 0; idet < kDet; idet++)
    {
      fHdetIso[idet] = NULL ;
      for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
	{
	  fHsmIso[idet][ismn] = NULL ;
	  for(Int_t jrow = 0; jrow < kMaxRow; jrow++)
	    {
	      for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
		{
		  fGainFact[idet][ismn][jrow][kcol] = 0.0;
		  fHadcIso[idet][ismn][jrow][kcol]  = NULL;
		}
	    }
	}
    }
}
// ------------------------------------------------------------------------ //
AliPMDCalibrator::AliPMDCalibrator(const AliPMDCalibrator &pmdcalibrator):
  fCalibGain(new AliPMDCalibData()), fCalibPed(new AliPMDPedestal())
{
  for(Int_t idet = 0; idet < 2; idet++)
    {
      fHdetIso[idet] = pmdcalibrator.fHdetIso[idet] ;
      for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
	{
	  fHsmIso[idet][ismn] = pmdcalibrator.fHsmIso[idet][ismn] ;
	  for(Int_t jrow = 0; jrow < kMaxRow; jrow++)
	    {
	      for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
		{
		  fGainFact[idet][ismn][jrow][kcol] = pmdcalibrator.fGainFact[idet][ismn][jrow][kcol];
		  fHadcIso[idet][ismn][jrow][kcol]  = pmdcalibrator.fHadcIso[idet][ismn][jrow][kcol];
		}
	    }
	}
    }
  
}
// ------------------------------------------------------------------------ //
AliPMDCalibrator &AliPMDCalibrator::operator=(const AliPMDCalibrator &pmdcalibrator)
{
  if(this != &pmdcalibrator)
    {
      for(Int_t idet = 0; idet < kDet; idet++)
	{
	  fHdetIso[idet] = pmdcalibrator.fHdetIso[idet] ;
	  for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
	    {
	      fHsmIso[idet][ismn] = pmdcalibrator.fHsmIso[idet][ismn] ;
	      for(Int_t jrow = 0; jrow < kMaxRow;jrow++)
		{
		  for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
		    {
		      fGainFact[idet][ismn][jrow][kcol] =
			pmdcalibrator.fGainFact[idet][ismn][jrow][kcol];
		      fHadcIso[idet][ismn][jrow][kcol]  =
			pmdcalibrator.fHadcIso[idet][ismn][jrow][kcol];
		    }
		}
	    }
	}
    }
  return *this;
}
// ------------------------------------------------------------------------ //
AliPMDCalibrator::~AliPMDCalibrator()
{
  // destructor
  for (Int_t i=0; i<kDet; i++)
    delete fHdetIso[i] ;
  for (Int_t i=0; i<kDet; i++)
    for (Int_t j=0; j<kMaxSMN; j++)
      delete fHsmIso[i][j] ;
  for (Int_t i=0; i<kDet; i++)
    for (Int_t j=0; j<kMaxSMN; j++)
      for (Int_t k=0; k<kMaxRow; k++)
	for (Int_t l=0; l<kMaxCol; l++)
	  delete fHadcIso[i][j][k][l] ;
  delete fCalibGain;
  delete fCalibPed;
}
// ------------------------------------------------------------------------ //

void AliPMDCalibrator::Exec()
{
  // reads parameters and does the calibration
  CalculateIsoCell() ;
}
// ------------------------------------------------------------------------ //

void AliPMDCalibrator::Init()
{
  // intializes everything
  char hname[50];
  char hname2[50];
  char hname24[50];
  char hnameiso[120];
  char htitle1[120];
  char htitle2[120];
  
  for(Int_t d = 0; d < kDet; d++) {
    snprintf(hname2,50,"Isolated cell adc for Det Plane %d",d);
    fHdetIso[d]= new TH1F(hname2,htitle2,100,0,1000);
    for(Int_t i1 = 0; i1 < kMaxSMN; i1++) {
      snprintf(hname,50,"det_%d_iso_sm_%2d",d,i1);
      snprintf(hname24,50,"det_%d_iso_sm_%2d",d,i1);
      fHsmIso[d][i1]= new TH1F(hname,hname24,100,0,1000);
      for(Int_t j1 = 0; j1 < kMaxRow; j1++) {
	for(Int_t k1 = 0; k1 < kMaxCol; k1++) {
	  snprintf(hnameiso,120,"Isolated Cell ADC for det_%d_cell_sm%d_row%d_col%d"
		  ,d,i1,j1,k1);
	  snprintf(htitle1,120,"Isolated Cell ADC for det_%d_cell_sm%d_row%d_col%d"
		  ,d,i1,j1,k1);
	  
	  TObject *old=gDirectory->GetList()->FindObject(hnameiso);
	  if (old) gDirectory->GetList()->Remove(old);
	  fHadcIso[d][i1][j1][k1] = new TH1F(hnameiso,htitle1,100,0.,1000.);
	}
      }
    }
  }
  
}

// ------------------------------------------------------------------------ //

void AliPMDCalibrator::CalculateIsoCell()
{
  // Calculates the ADC of isolated cell

  TObjArray pmdddlcont;

  const Int_t kCellNeighbour = 6;

  Int_t neibx[6] = {1,0,-1,-1,0,1};
  Int_t neiby[6] = {0,1,1,0,-1,-1};

  Int_t id1 = 0,jd1 = 0;        //neighbour row/col
  Int_t countisocell = 0 ;      //number of isilated cell
  Int_t isocount = 0;           //number of neighbours with 0 signal
  Int_t maxhit = 0;
  Int_t nhit[kDet][kMaxSMN];
  Int_t d1[kDet][kMaxSMN][kMaxRow][kMaxCol];
  Int_t nhitcell[kDet][kMaxSMN][kMaxRow][kMaxCol];
  
  for(Int_t idet = 0; idet < kDet; idet++)
    {
      for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
        {
	  nhit[idet][ismn] = 0;
	  for(Int_t irow = 0; irow < kMaxRow; irow++)
	    {
	      for(Int_t icol = 0; icol < kMaxCol; icol++)
		{
		  d1[idet][ismn][irow][icol] = 0;
		  nhitcell[idet][ismn][irow][icol] = 0;
		}
	    }
	}
    }

  Float_t tempnhit1[kDet][kMaxSMN];
  Float_t tempnhit2[kDet][kMaxSMN];
  Float_t meannhit[kDet][kMaxSMN];
  Float_t meanSqnhit[kDet][kMaxSMN];
  Float_t sigmanhit[kDet][kMaxSMN];
  Float_t count[kDet][kMaxSMN];
  Float_t nhitcut[kDet][kMaxSMN];

  for (Int_t idet = 0; idet < kDet; idet++)
  {
      for (Int_t ismn = 0; ismn < kMaxSMN; ismn++)
      {
	  tempnhit1[idet][ismn]  = 0.;
	  tempnhit2[idet][ismn]  = 0.;
	  meannhit[idet][ismn]   = 0.;
	  meanSqnhit[idet][ismn] = 0.;
	  sigmanhit[idet][ismn]  = 0.;
	  count[idet][ismn]      = 0.;
	  nhitcut[idet][ismn]    = 0.;
      }
  }


  //accessing raw data
  AliRawReaderFile reader(".");
  AliPMDRawStream stream(&reader);
  while(reader.NextEvent())
    { 
      // New PMD Reader is plugged in
	Int_t iddl = -1;
	while ((iddl = stream.DdlData(&pmdddlcont)) >=0) {

	    //reader.Select("PMD", iddl, iddl);
	    //stream.DdlData(iddl,&pmdddlcont);
	  Int_t ientries = pmdddlcont.GetEntries();
	  for (Int_t ient = 0; ient < ientries; ient++)
	    {
	      AliPMDddldata *pmdddl = 
		(AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
	      Int_t idet = pmdddl->GetDetector();
	      Int_t ismn = pmdddl->GetSMN();
	      Int_t irow = pmdddl->GetRow();
	      Int_t icol = pmdddl->GetColumn();
	      Float_t isig1 = pmdddl->GetSignal();
	      // Pedestal Subtraction
	      Int_t   pedmeanrms = 
		fCalibPed->GetPedMeanRms(idet,ismn,irow,icol);
	      Int_t   pedrms1    = (Int_t) pedmeanrms%100;
	      Float_t pedrms     = (Float_t)pedrms1/10.;
	      Float_t pedmean    = (Float_t) (pedmeanrms - pedrms1)/1000.0;
	      Float_t isig = isig1 - (pedmean + 3.0*pedrms);
	      if (isig>0)
		{
		  d1[idet][ismn][irow][icol] = (Int_t)isig;
		  nhitcell[idet][ismn][irow][icol] += 1; 
		}
	    }//ient loop
	  pmdddlcont.Clear();
	}//iddl loop
      maxhit = 0;
      for(Int_t idet=0; idet < kDet; idet++)
	{
	  for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
	    {
	      for(Int_t irow = 0; irow < kMaxRow; irow++)
		{
		  for(Int_t icol = 0; icol < kMaxCol; icol++)
		    {
		      if(d1[idet][ismn][irow][icol] > 0)
			{
			  isocount = 0;
			  for(Int_t ii = 0; ii < kCellNeighbour; ii++)
			    {
			      id1 = irow + neibx[ii];
			      jd1 = icol + neiby[ii];
			      if(d1[idet][ismn][id1][jd1] == 0)
				{
				  isocount++;
				  if(isocount == kCellNeighbour)
				    {
				      countisocell++;
				      maxhit++;
                                      fHdetIso[idet]->
					Fill(d1[idet][ismn][irow][icol]);
				      fHsmIso[idet][ismn]->
					Fill(d1[idet][ismn][irow][icol]);
				      fHadcIso[idet][ismn][irow][icol]->
					Fill(d1[idet][ismn][irow][icol]);
				    }
				}
			    }  // neigh cell cond.
			}     // d>0 cond.
		    }
		}
	    }
	} //det
    }//event loop

  // Mean and Sigma Calculations
  for(Int_t idet=0; idet < kDet; idet++){
      for(Int_t ismn = 0; ismn < kMaxSMN; ismn++){
	  for(Int_t irow = 0; irow < kMaxRow; irow++){
	      for(Int_t icol = 0; icol < kMaxCol; icol++){
		  if(nhitcell[idet][ismn][irow][icol]>0){
		      count[idet][ismn] += 1;
		      tempnhit1[idet][ismn] += nhitcell[idet][ismn][irow][icol];
		      tempnhit2[idet][ismn] += nhitcell[idet][ismn][irow][icol]
			  *nhitcell[idet][ismn][irow][icol];
		  }
	      }
	  }
      }
  }//det loop

  //cout<<"nhit cell = "<<idet<<"  "<<ismn<<"  "
  //  <<irow<<"  "<<icol<<"  "<<nhitcell[idet][ismn][irow][icol]<<endl;
  //count[idet][ismn] += 1;

  for(Int_t i=0; i < kDet; i++)
    {
      for(Int_t j=0; j < kMaxSMN; j++)
	{
	  if(count[i][j] > 0.)
	    {
		meannhit[i][j]   = tempnhit1[i][j]/count[i][j];
		meanSqnhit[i][j] = tempnhit2[i][j]/count[i][j]; 
		sigmanhit[i][j]  = sqrt(meanSqnhit[i][j]-
				      (meannhit[i][j]*meannhit[i][j]));
		nhitcut[i][j]    = 3*sigmanhit[i][j] + meannhit[i][j];
	    }
	}
    }


  Double_t histdetMean[kDet];
  Double_t histMean[kDet][kMaxSMN];
  Double_t isoMean[kDet][kMaxSMN][kMaxRow][kMaxCol];
  Double_t smNormFactor[kDet][kMaxSMN];
  
  for(Int_t det1 = 0; det1 < kDet; det1++)
    {
      histdetMean[det1]= fHdetIso[det1]->GetMean();
      for(Int_t i1 = 0; i1 < kMaxSMN; i1++)
	{
	  histMean[det1][i1]= fHsmIso[det1][i1]->GetMean();
	  if(histMean[det1][i1]>0.0 && histdetMean[det1]>0.0)
	    {
	      smNormFactor[det1][i1]= histdetMean[det1]/histMean[det1][i1];
	    }
	  for(Int_t j1 = 0; j1 < kMaxRow; j1++)
	    {
	      for(Int_t k1 = 0; k1 < kMaxCol; k1++)
		{
		 if(nhitcell[det1][i1][j1][k1]> nhitcut[det1][i1]) fGainFact[det1][i1][j1][k1]=-99.0;
		   if(nhitcell[det1][i1][j1][k1]< nhitcut[det1][i1])
		     {
		      isoMean[det1][i1][j1][k1]=fHadcIso[det1][i1][j1][k1]->
			GetMean();
		      if(isoMean[det1][i1][j1][k1]>0.0 && histMean[det1][i1]>0.0)
			{
			  fGainFact[det1][i1][j1][k1]=
			    isoMean[det1][i1][j1][k1]/(histMean[det1][i1]*
						     smNormFactor[det1][i1]);

			}
		    }   
 		   Float_t gain = fGainFact[det1][i1][j1][k1];
 		   fCalibGain->SetGainFact(det1,i1,j1,k1,gain);
		}
                          
 			
	    }
	}
    }
}//CalculateIsoCell()

// ------------------------------------------------------------------------ //

Bool_t AliPMDCalibrator::Store()
{
  AliCDBManager *man = AliCDBManager::Instance();
  //man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
  if(!man->IsDefaultStorageSet()) return kFALSE;
  AliCDBId id("PMD/Calib/Gain",0,999999999);
  AliCDBMetaData md;
  md.SetBeamPeriod(0);
  md.SetComment("Test");
  
  printf("\n\n\n fCalibData\n");
  //fCalibData->Print(0);
  //printf("\n\n\n fCalibData\n");
  
  Bool_t result = man->Put(fCalibGain,id,&md);

  return result;
}
Commit	Line	Data
06810de6	1	//////////////////////////////////////////////////////////////////////////////
	2	//
	3	// * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	4	// * *
	5	// * Author: The ALICE Off-line Project. *
	6	// * Contributors are mentioned in the code where appropriate. *
	7	// * *
	8	// * Permission to use, copy, modify and distribute this software and its *
	9	// * documentation strictly for non-commercial purposes is hereby granted *
	10	// * without fee, provided that the above copyright notice appears in all *
	11	// * copies and that both the copyright notice and this permission notice *
	12	// * appear in the supporting documentation. The authors make no claims *
	13	// * about the suitability of this software for any purpose. It is *
	14	// * provided "as is" without express or implied warranty. *
	15	// **************************************************************************/
	16	//
	17	//////////////////////////////////////////////////////////////////////////////
06810de6	18
	19	#include "TF1.h"
	20	#include "TFile.h"
	21	#include "TObjString.h"
	22	#include "TROOT.h"
	23	#include "TClonesArray.h"
	24	#include "TH1F.h"
df42ab21	25	#include "TObjArray.h"
54a2e959	26	#include "TMath.h"
06810de6	27
	28	// --- Standard library ---
	29
	30	// --- AliRoot header files ---
	31	#include "AliLog.h"
df42ab21	32	#include "AliRawReaderFile.h"
06810de6	33	#include "AliPMDCalibrator.h"
	34	#include "AliRawReaderDate.h"
	35	#include "AliPMDRawStream.h"
06810de6	36	#include "AliPMDCalibData.h"
df42ab21	37	#include "AliPMDddldata.h"
06810de6	38	#include "AliCDBManager.h"
	39	#include "AliCDBId.h"
	40	#include "AliCDBMetaData.h"
df42ab21	41	#include "AliDAQ.h"
06810de6	42
54a2e959	43	//#include "AliPMDCleanNoise.h"
	44	//#include "AliPMDCleaner.h"
	45
	46	#include "AliPMDPedestal.h"
	47	#include "AliCDBManager.h"
	48	#include "AliCDBEntry.h"
	49
06810de6	50	ClassImp(AliPMDCalibrator)
06810de6	51
06810de6	52
a48edddd	53	AliPMDCalibrator::AliPMDCalibrator():
54a2e959	54	fCalibGain(new AliPMDCalibData()), fCalibPed(new AliPMDPedestal())
06810de6	55	{
06810de6	56	// Standard Constructor
35535af7	57	for(Int_t idet = 0; idet < kDet; idet++)
06810de6	58	{
54a2e959	59	fHdetIso[idet] = NULL ;
35535af7	60	for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
06810de6	61	{
35535af7	62	fHsmIso[idet][ismn] = NULL ;
35535af7	63	for(Int_t jrow = 0; jrow < kMaxRow; jrow++)
06810de6	64	{
35535af7	65	for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
06810de6	66	{
35535af7	67	fGainFact[idet][ismn][jrow][kcol] = 0.0;
35535af7	68	fHadcIso[idet][ismn][jrow][kcol] = NULL;
06810de6	69	}
	70	}
	71	}
	72	}
	73	}
	74	// ------------------------------------------------------------------------ //
a48edddd	75	AliPMDCalibrator::AliPMDCalibrator(const AliPMDCalibrator &pmdcalibrator):
54a2e959	76	fCalibGain(new AliPMDCalibData()), fCalibPed(new AliPMDPedestal())
a48edddd	77	{
35535af7	78	for(Int_t idet = 0; idet < 2; idet++)
a48edddd	79	{
54a2e959	80	fHdetIso[idet] = pmdcalibrator.fHdetIso[idet] ;
35535af7	81	for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
a48edddd	82	{
35535af7	83	fHsmIso[idet][ismn] = pmdcalibrator.fHsmIso[idet][ismn] ;
35535af7	84	for(Int_t jrow = 0; jrow < kMaxRow; jrow++)
a48edddd	85	{
35535af7	86	for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
a48edddd	87	{
35535af7	88	fGainFact[idet][ismn][jrow][kcol] = pmdcalibrator.fGainFact[idet][ismn][jrow][kcol];
35535af7	89	fHadcIso[idet][ismn][jrow][kcol] = pmdcalibrator.fHadcIso[idet][ismn][jrow][kcol];
a48edddd	90	}
	91	}
	92	}
	93	}
54a2e959	94
a48edddd	95	}
	96	// ------------------------------------------------------------------------ //
	97	AliPMDCalibrator &AliPMDCalibrator::operator=(const AliPMDCalibrator &pmdcalibrator)
	98	{
	99	if(this != &pmdcalibrator)
	100	{
35535af7	101	for(Int_t idet = 0; idet < kDet; idet++)
a48edddd	102	{
54a2e959	103	fHdetIso[idet] = pmdcalibrator.fHdetIso[idet] ;
35535af7	104	for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
a48edddd	105	{
35535af7	106	fHsmIso[idet][ismn] = pmdcalibrator.fHsmIso[idet][ismn] ;
35535af7	107	for(Int_t jrow = 0; jrow < kMaxRow;jrow++)
a48edddd	108	{
35535af7	109	for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
a48edddd	110	{
35535af7	111	fGainFact[idet][ismn][jrow][kcol] =
	112	pmdcalibrator.fGainFact[idet][ismn][jrow][kcol];
	113	fHadcIso[idet][ismn][jrow][kcol] =
	114	pmdcalibrator.fHadcIso[idet][ismn][jrow][kcol];
a48edddd	115	}
	116	}
	117	}
	118	}
	119	}
	120	return *this;
	121	}
	122	// ------------------------------------------------------------------------ //
06810de6	123	AliPMDCalibrator::~AliPMDCalibrator()
06810de6	124	{
54a2e959	125	// destructor
63424920	126	for (Int_t i=0; i<kDet; i++)
	127	delete fHdetIso[i] ;
	128	for (Int_t i=0; i<kDet; i++)
	129	for (Int_t j=0; j<kMaxSMN; j++)
	130	delete fHsmIso[i][j] ;
	131	for (Int_t i=0; i<kDet; i++)
	132	for (Int_t j=0; j<kMaxSMN; j++)
	133	for (Int_t k=0; k<kMaxRow; k++)
	134	for (Int_t l=0; l<kMaxCol; l++)
	135	delete fHadcIso[i][j][k][l] ;
35535af7	136	delete fCalibGain;
54a2e959	137	delete fCalibPed;
06810de6	138	}
	139	// ------------------------------------------------------------------------ //
	140
	141	void AliPMDCalibrator::Exec()
	142	{
	143	// reads parameters and does the calibration
	144	CalculateIsoCell() ;
06810de6	145	}
	146	// ------------------------------------------------------------------------ //
	147
	148	void AliPMDCalibrator::Init()
	149	{
	150	// intializes everything
4e2efa25	151	char hname[50];
	152	char hname2[50];
	153	char hname24[50];
06810de6	154	char hnameiso[120];
06810de6	155	char htitle1[120];
4e2efa25	156	char htitle2[120];
54a2e959	157
35535af7	158	for(Int_t d = 0; d < kDet; d++) {
4e2efa25	159	snprintf(hname2,50,"Isolated cell adc for Det Plane %d",d);
54a2e959	160	fHdetIso[d]= new TH1F(hname2,htitle2,100,0,1000);
35535af7	161	for(Int_t i1 = 0; i1 < kMaxSMN; i1++) {
4e2efa25	162	snprintf(hname,50,"det_%d_iso_sm_%2d",d,i1);
4e2efa25	163	snprintf(hname24,50,"det_%d_iso_sm_%2d",d,i1);
06810de6	164	fHsmIso[d][i1]= new TH1F(hname,hname24,100,0,1000);
	165	for(Int_t j1 = 0; j1 < kMaxRow; j1++) {
	166	for(Int_t k1 = 0; k1 < kMaxCol; k1++) {
4e2efa25	167	snprintf(hnameiso,120,"Isolated Cell ADC for det_%d_cell_sm%d_row%d_col%d"
06810de6	168	,d,i1,j1,k1);
4e2efa25	169	snprintf(htitle1,120,"Isolated Cell ADC for det_%d_cell_sm%d_row%d_col%d"
06810de6	170	,d,i1,j1,k1);
	171
	172	TObject *old=gDirectory->GetList()->FindObject(hnameiso);
	173	if (old) gDirectory->GetList()->Remove(old);
54a2e959	174	fHadcIso[d][i1][j1][k1] = new TH1F(hnameiso,htitle1,100,0.,1000.);
06810de6	175	}
	176	}
	177	}
	178	}
	179
	180	}
	181
	182	// ------------------------------------------------------------------------ //
	183
	184	void AliPMDCalibrator::CalculateIsoCell()
	185	{
	186	// Calculates the ADC of isolated cell
	187
df42ab21	188	TObjArray pmdddlcont;
722ccc67	189
06810de6	190	const Int_t kCellNeighbour = 6;
	191
	192	Int_t neibx[6] = {1,0,-1,-1,0,1};
	193	Int_t neiby[6] = {0,1,1,0,-1,-1};
54a2e959	194
be8b7039	195	Int_t id1 = 0,jd1 = 0; //neighbour row/col
	196	Int_t countisocell = 0 ; //number of isilated cell
	197	Int_t isocount = 0; //number of neighbours with 0 signal
	198	Int_t maxhit = 0;
54a2e959	199	Int_t nhit[kDet][kMaxSMN];
be8b7039	200	Int_t d1[kDet][kMaxSMN][kMaxRow][kMaxCol];
54a2e959	201	Int_t nhitcell[kDet][kMaxSMN][kMaxRow][kMaxCol];
06810de6	202
06810de6	203	for(Int_t idet = 0; idet < kDet; idet++)
	204	{
	205	for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
	206	{
54a2e959	207	nhit[idet][ismn] = 0;
	208	for(Int_t irow = 0; irow < kMaxRow; irow++)
	209	{
	210	for(Int_t icol = 0; icol < kMaxCol; icol++)
	211	{
	212	d1[idet][ismn][irow][icol] = 0;
	213	nhitcell[idet][ismn][irow][icol] = 0;
	214	}
06810de6	215	}
	216	}
	217	}
54a2e959	218
	219	Float_t tempnhit1[kDet][kMaxSMN];
	220	Float_t tempnhit2[kDet][kMaxSMN];
	221	Float_t meannhit[kDet][kMaxSMN];
	222	Float_t meanSqnhit[kDet][kMaxSMN];
	223	Float_t sigmanhit[kDet][kMaxSMN];
	224	Float_t count[kDet][kMaxSMN];
	225	Float_t nhitcut[kDet][kMaxSMN];
	226
	227	for (Int_t idet = 0; idet < kDet; idet++)
	228	{
	229	for (Int_t ismn = 0; ismn < kMaxSMN; ismn++)
	230	{
	231	tempnhit1[idet][ismn] = 0.;
	232	tempnhit2[idet][ismn] = 0.;
	233	meannhit[idet][ismn] = 0.;
	234	meanSqnhit[idet][ismn] = 0.;
	235	sigmanhit[idet][ismn] = 0.;
	236	count[idet][ismn] = 0.;
	237	nhitcut[idet][ismn] = 0.;
	238	}
	239	}
	240
	241
06810de6	242	//accessing raw data
	243	AliRawReaderFile reader(".");
	244	AliPMDRawStream stream(&reader);
	245	while(reader.NextEvent())
54a2e959	246	{
df42ab21	247	// New PMD Reader is plugged in
722ccc67	248	Int_t iddl = -1;
	249	while ((iddl = stream.DdlData(&pmdddlcont)) >=0) {
	250
	251	//reader.Select("PMD", iddl, iddl);
	252	//stream.DdlData(iddl,&pmdddlcont);
df42ab21	253	Int_t ientries = pmdddlcont.GetEntries();
	254	for (Int_t ient = 0; ient < ientries; ient++)
	255	{
54a2e959	256	AliPMDddldata *pmdddl =
54a2e959	257	(AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
df42ab21	258	Int_t idet = pmdddl->GetDetector();
df42ab21	259	Int_t ismn = pmdddl->GetSMN();
df42ab21	260	Int_t irow = pmdddl->GetRow();
df42ab21	261	Int_t icol = pmdddl->GetColumn();
54a2e959	262	Float_t isig1 = pmdddl->GetSignal();
	263	// Pedestal Subtraction
	264	Int_t pedmeanrms =
	265	fCalibPed->GetPedMeanRms(idet,ismn,irow,icol);
37edc588	266	Int_t pedrms1 = (Int_t) pedmeanrms%100;
54a2e959	267	Float_t pedrms = (Float_t)pedrms1/10.;
	268	Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
	269	Float_t isig = isig1 - (pedmean + 3.0*pedrms);
df42ab21	270	if (isig>0)
df42ab21	271	{
54a2e959	272	d1[idet][ismn][irow][icol] = (Int_t)isig;
121747cb	273	nhitcell[idet][ismn][irow][icol] += 1;
df42ab21	274	}
54a2e959	275	}//ient loop
df42ab21	276	pmdddlcont.Clear();
54a2e959	277	}//iddl loop
06810de6	278	maxhit = 0;
06810de6	279	for(Int_t idet=0; idet < kDet; idet++)
	280	{
	281	for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
	282	{
	283	for(Int_t irow = 0; irow < kMaxRow; irow++)
	284	{
	285	for(Int_t icol = 0; icol < kMaxCol; icol++)
	286	{
06810de6	287	if(d1[idet][ismn][irow][icol] > 0)
	288	{
	289	isocount = 0;
	290	for(Int_t ii = 0; ii < kCellNeighbour; ii++)
	291	{
	292	id1 = irow + neibx[ii];
	293	jd1 = icol + neiby[ii];
	294	if(d1[idet][ismn][id1][jd1] == 0)
	295	{
	296	isocount++;
	297	if(isocount == kCellNeighbour)
	298	{
	299	countisocell++;
06810de6	300	maxhit++;
54a2e959	301	fHdetIso[idet]->
	302	Fill(d1[idet][ismn][irow][icol]);
	303	fHsmIso[idet][ismn]->
	304	Fill(d1[idet][ismn][irow][icol]);
	305	fHadcIso[idet][ismn][irow][icol]->
	306	Fill(d1[idet][ismn][irow][icol]);
06810de6	307	}
	308	}
	309	} // neigh cell cond.
	310	} // d>0 cond.
	311	}
	312	}
	313	}
54a2e959	314	} //det
	315	}//event loop
	316
	317	// Mean and Sigma Calculations
	318	for(Int_t idet=0; idet < kDet; idet++){
	319	for(Int_t ismn = 0; ismn < kMaxSMN; ismn++){
	320	for(Int_t irow = 0; irow < kMaxRow; irow++){
	321	for(Int_t icol = 0; icol < kMaxCol; icol++){
	322	if(nhitcell[idet][ismn][irow][icol]>0){
	323	count[idet][ismn] += 1;
	324	tempnhit1[idet][ismn] += nhitcell[idet][ismn][irow][icol];
	325	tempnhit2[idet][ismn] += nhitcell[idet][ismn][irow][icol]
	326	*nhitcell[idet][ismn][irow][icol];
	327	}
	328	}
	329	}
	330	}
	331	}//det loop
	332
	333	//cout<<"nhit cell = "<<idet<<" "<<ismn<<" "
	334	// <<irow<<" "<<icol<<" "<<nhitcell[idet][ismn][irow][icol]<<endl;
	335	//count[idet][ismn] += 1;
	336
	337	for(Int_t i=0; i < kDet; i++)
	338	{
	339	for(Int_t j=0; j < kMaxSMN; j++)
	340	{
	341	if(count[i][j] > 0.)
	342	{
	343	meannhit[i][j] = tempnhit1[i][j]/count[i][j];
	344	meanSqnhit[i][j] = tempnhit2[i][j]/count[i][j];
	345	sigmanhit[i][j] = sqrt(meanSqnhit[i][j]-
	346	(meannhit[i][j]*meannhit[i][j]));
	347	nhitcut[i][j] = 3*sigmanhit[i][j] + meannhit[i][j];
	348	}
	349	}
06810de6	350	}
54a2e959	351
	352
	353	Double_t histdetMean[kDet];
35535af7	354	Double_t histMean[kDet][kMaxSMN];
35535af7	355	Double_t isoMean[kDet][kMaxSMN][kMaxRow][kMaxCol];
54a2e959	356	Double_t smNormFactor[kDet][kMaxSMN];
54a2e959	357
78fc1b96	358	for(Int_t det1 = 0; det1 < kDet; det1++)
06810de6	359	{
78fc1b96	360	histdetMean[det1]= fHdetIso[det1]->GetMean();
35535af7	361	for(Int_t i1 = 0; i1 < kMaxSMN; i1++)
06810de6	362	{
78fc1b96	363	histMean[det1][i1]= fHsmIso[det1][i1]->GetMean();
78fc1b96	364	if(histMean[det1][i1]>0.0 && histdetMean[det1]>0.0)
54a2e959	365	{
78fc1b96	366	smNormFactor[det1][i1]= histdetMean[det1]/histMean[det1][i1];
54a2e959	367	}
35535af7	368	for(Int_t j1 = 0; j1 < kMaxRow; j1++)
06810de6	369	{
35535af7	370	for(Int_t k1 = 0; k1 < kMaxCol; k1++)
06810de6	371	{
121747cb	372	if(nhitcell[det1][i1][j1][k1]> nhitcut[det1][i1]) fGainFact[det1][i1][j1][k1]=-99.0;
	373	if(nhitcell[det1][i1][j1][k1]< nhitcut[det1][i1])
	374	{
78fc1b96	375	isoMean[det1][i1][j1][k1]=fHadcIso[det1][i1][j1][k1]->
54a2e959	376	GetMean();
78fc1b96	377	if(isoMean[det1][i1][j1][k1]>0.0 && histMean[det1][i1]>0.0)
54a2e959	378	{
78fc1b96	379	fGainFact[det1][i1][j1][k1]=
	380	isoMean[det1][i1][j1][k1]/(histMean[det1][i1]*
	381	smNormFactor[det1][i1]);
121747cb	382
54a2e959	383	}
121747cb	384	}
	385	Float_t gain = fGainFact[det1][i1][j1][k1];
	386	fCalibGain->SetGainFact(det1,i1,j1,k1,gain);
06810de6	387	}
121747cb	388
121747cb	389
06810de6	390	}
	391	}
	392	}
54a2e959	393	}//CalculateIsoCell()
54a2e959	394
06810de6	395	// ------------------------------------------------------------------------ //
54a2e959	396
06810de6	397	Bool_t AliPMDCalibrator::Store()
	398	{
	399	AliCDBManager *man = AliCDBManager::Instance();
162637e4	400	//man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
0dd3d6f9	401	if(!man->IsDefaultStorageSet()) return kFALSE;
121747cb	402	AliCDBId id("PMD/Calib/Gain",0,999999999);
06810de6	403	AliCDBMetaData md;
06810de6	404	md.SetBeamPeriod(0);
06810de6	405	md.SetComment("Test");
	406
	407	printf("\n\n\n fCalibData\n");
	408	//fCalibData->Print(0);
	409	//printf("\n\n\n fCalibData\n");
	410
35535af7	411	Bool_t result = man->Put(fCalibGain,id,&md);
06810de6	412
	413	return result;
	414	}
	415