[u/mrichter/AliRoot.git] / TRD / Cal / AliTRDCalROC.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$ */

///////////////////////////////////////////////////////////////////////////////
//                                                                           //
//  Calibration base class for a single ROC                                  //
//  Contains one UShort_t value per pad                                      //
//  However, values are set and get as float, there are stored internally as //
//  (UShort_t) value * 10000                                                 //
//                                                                           //
///////////////////////////////////////////////////////////////////////////////

#include <TStyle.h>
#include <TMath.h>
#include <TH2F.h>
#include <TH1F.h>

#include "AliMathBase.h"

#include "AliTRDCalROC.h"

ClassImp(AliTRDCalROC)

//_____________________________________________________________________________
AliTRDCalROC::AliTRDCalROC()
  :TObject()
  ,fPla(0)
  ,fCha(0)
  ,fNrows(0)
  ,fNcols(0)
  ,fNchannels(0)
  ,fData(0)
{
  //
  // Default constructor
  //

}

//_____________________________________________________________________________
AliTRDCalROC::AliTRDCalROC(Int_t p, Int_t c)
  :TObject()
  ,fPla(p)
  ,fCha(c)
  ,fNrows(0)
  ,fNcols(144)
  ,fNchannels(0)
  ,fData(0)
{
  //
  // Constructor that initializes a given pad plane type
  //

  //
  // The pad plane parameter
  //
  switch (p) {
  case 0:
    if (c == 2) {
      // L0C0 type
      fNrows        =  12;
    }
    else {
      // L0C1 type
      fNrows        =  16;
    }
    break;
  case 1:
    if (c == 2) {
      // L1C0 type
      fNrows        =  12;
    }
    else {
      // L1C1 type
      fNrows        =  16;
    }
    break;
  case 2:
    if (c == 2) {
      // L2C0 type
      fNrows        =  12;
    }
    else {
      // L2C1 type
      fNrows        =  16;
    }
    break;
  case 3:
    if (c == 2) {
      // L3C0 type
      fNrows        =  12;
    }
    else {
      // L3C1 type
      fNrows        =  16;
    }
    break;
  case 4:
    if (c == 2) {
      // L4C0 type
      fNrows        =  12;
    }
    else {
      // L4C1 type
      fNrows        =  16;
    }
    break;
  case 5:
    if (c == 2) {
      // L5C0 type
      fNrows        =  12;
    }
    else {
      // L5C1 type
      fNrows        =  16;
    }
    break;
  };

  fNchannels = fNrows * fNcols;
  if (fNchannels != 0) {
    fData = new UShort_t[fNchannels];
  }

  for (Int_t i=0; i<fNchannels; ++i) {
    fData[i] = 0;
  }

}

//_____________________________________________________________________________
AliTRDCalROC::AliTRDCalROC(const AliTRDCalROC &c)
  :TObject(c)
  ,fPla(c.fPla)
  ,fCha(c.fCha)
  ,fNrows(c.fNrows)
  ,fNcols(c.fNcols)
  ,fNchannels(c.fNchannels)
  ,fData(0)
{
  //
  // AliTRDCalROC copy constructor
  //

  Int_t iBin = 0;

  fData = new UShort_t[fNchannels];
  for (iBin = 0; iBin < fNchannels; iBin++) {
    fData[iBin] = ((AliTRDCalROC &) c).fData[iBin];
  }

}

//_____________________________________________________________________________
AliTRDCalROC::~AliTRDCalROC()
{
  //
  // AliTRDCalROC destructor
  //

  if (fData) {
    delete [] fData;
    fData = 0;
  }

}

//_____________________________________________________________________________
AliTRDCalROC &AliTRDCalROC::operator=(const AliTRDCalROC &c)
{
  //
  // Assignment operator
  //

  if (this != &c) ((AliTRDCalROC &) c).Copy(*this);
  return *this;

}

//_____________________________________________________________________________
void AliTRDCalROC::Copy(TObject &c) const
{
  //
  // Copy function
  //

  ((AliTRDCalROC &) c).fPla          = fPla;
  ((AliTRDCalROC &) c).fCha          = fCha;

  ((AliTRDCalROC &) c).fNrows        = fNrows;
  ((AliTRDCalROC &) c).fNcols        = fNcols;

  Int_t iBin = 0;

  ((AliTRDCalROC &) c).fNchannels = fNchannels;

  if (((AliTRDCalROC &) c).fData) delete [] ((AliTRDCalROC &) c).fData;
  ((AliTRDCalROC &) c).fData = new UShort_t[fNchannels];
  for (iBin = 0; iBin < fNchannels; iBin++) {
    ((AliTRDCalROC &) c).fData[iBin] = fData[iBin];
  }
  
  TObject::Copy(c);

}

//___________________________________________________________________________________
Double_t AliTRDCalROC::GetMean(AliTRDCalROC* const outlierROC) const
{
  //
  // Calculate the mean
  //

   Double_t *ddata = new Double_t[fNchannels];
   Int_t nPoints = 0;
   for (Int_t i=0;i<fNchannels;i++) {
      if ((!outlierROC) || (!(outlierROC->GetValue(i)))) {
	if(fData[i] > 0.000000000000001){
         ddata[nPoints]= (Double_t) fData[i]/10000;
         nPoints++;
	}
      }
   }
   Double_t mean = TMath::Mean(nPoints,ddata);
   delete [] ddata;
   return mean;
}

//_______________________________________________________________________________________
Double_t AliTRDCalROC::GetMedian(AliTRDCalROC* const outlierROC) const
{
  //
  // Calculate the median
  //

  Double_t *ddata = new Double_t[fNchannels];
   Int_t nPoints = 0;
   for (Int_t i=0;i<fNchannels;i++) {
       if ((!outlierROC) || (!(outlierROC->GetValue(i)))) {
	 if(fData[i] > 0.000000000000001){         
	   ddata[nPoints]= (Double_t) fData[i]/10000;
	   nPoints++;
	 }
       }
   }
   Double_t mean = TMath::Median(nPoints,ddata);
   delete [] ddata;
   return mean;
}

//____________________________________________________________________________________________
Double_t AliTRDCalROC::GetRMS(AliTRDCalROC* const outlierROC) const
{
  //
  // Calculate the RMS
  //

  Double_t *ddata = new Double_t[fNchannels];
  Int_t nPoints = 0;
  for (Int_t i=0;i<fNchannels;i++) {
    if ((!outlierROC) || (!(outlierROC->GetValue(i)))) {
       if(fData[i] > 0.000000000000001){
         ddata[nPoints]= (Double_t)fData[i]/10000;
         nPoints++;
       }
    }
  }
  Double_t mean = TMath::RMS(nPoints,ddata);
  delete [] ddata;
  return mean;
}

//______________________________________________________________________________________________
Double_t AliTRDCalROC::GetLTM(Double_t *sigma, Double_t fraction, AliTRDCalROC* const outlierROC)
{
  //
  //  Calculate LTM mean and sigma
  //

  Double_t *ddata = new Double_t[fNchannels];
  Double_t mean=0, lsigma=0;
  UInt_t nPoints = 0;
  for (Int_t i=0;i<fNchannels;i++) {
     if (!outlierROC || !(outlierROC->GetValue(i))) {
       if(fData[i] > 0.000000000000001){
	 ddata[nPoints]= (Double_t) fData[i]/10000;
	 nPoints++;
       }
     }
  }
  Int_t hh = TMath::Min(TMath::Nint(fraction *nPoints), Int_t(nPoints));
  AliMathBase::EvaluateUni(nPoints,ddata, mean, lsigma, hh);
  if (sigma) *sigma=lsigma;
  delete [] ddata;
  return mean;
}

//___________________________________________________________________________________
Bool_t AliTRDCalROC::Add(Float_t c1)
{
  //
  // add constant
  //

  Bool_t result = kTRUE;
  for (Int_t  idata = 0; idata< fNchannels; idata++) {
    if(((GetValue(idata)+c1) <= 6.5535) && ((GetValue(idata)+c1) >= 0.0)) SetValue(idata,GetValue(idata)+c1);
    else {
      SetValue(idata,0.0);
      result = kFALSE;
    }
  }
  return result;
}

//_______________________________________________________________________________________
Bool_t AliTRDCalROC::Multiply(Float_t c1)
{
  //
  // multiply constant
  //

  Bool_t result = kTRUE;
  if(c1 < 0) return kFALSE;
  for (Int_t  idata = 0; idata< fNchannels; idata++) {
    if((GetValue(idata)*c1) <= 6.5535) SetValue(idata,GetValue(idata)*c1);
    else {
      SetValue(idata,0.0);
      result = kFALSE;
    }
  }
  return result;
}

//____________________________________________________________________________________________
Bool_t AliTRDCalROC::Add(const AliTRDCalROC * roc, Double_t c1)
{
  //
  // add values 
  //

  Bool_t result = kTRUE;
  for (Int_t  idata = 0; idata< fNchannels; idata++){
    if(((GetValue(idata)+roc->GetValue(idata)*c1) <= 6.5535) && 
       ((GetValue(idata)+roc->GetValue(idata)*c1) >= 0.0)) 
          SetValue(idata,GetValue(idata)+roc->GetValue(idata)*c1);
    else {
      SetValue(idata,0.0);
      result = kFALSE;
    }
  }
  return result;
}

//____________________________________________________________________________________________
Bool_t AliTRDCalROC::Multiply(const AliTRDCalROC*  roc) 
{
  //
  // multiply values - per by pad
  //

  Bool_t result = kTRUE;
  for (Int_t  idata = 0; idata< fNchannels; idata++){
    if((GetValue(idata)*roc->GetValue(idata)) <= 6.5535) 
      SetValue(idata,GetValue(idata)*roc->GetValue(idata));
    else {
      SetValue(idata,0.0);
      result = kFALSE;
    }
  }
  return result;
}

//______________________________________________________________________________________________
Bool_t AliTRDCalROC::Divide(const AliTRDCalROC*  roc) 
{
  //
  // divide values 
  //

  Bool_t result = kTRUE;
  Float_t kEpsilon=0.00000000000000001;
  for (Int_t  idata = 0; idata< fNchannels; idata++){
    if (TMath::Abs(roc->GetValue(idata))>kEpsilon){
      if((GetValue(idata)/roc->GetValue(idata)) <= 6.5535) 
        SetValue(idata,GetValue(idata)/roc->GetValue(idata));
      else {
	SetValue(idata,0.0);
	result = kFALSE;
      }
    }
    else result = kFALSE;
  }
  return result;
}
//______________________________________________________________________________________________
Bool_t AliTRDCalROC::Unfold() 
{
  //
  // Compute the mean value per pad col
  // Divide with this value each pad col
  // This is for the noise study
  // Return kFALSE if one or more of the pad col was not normalised
  //

  Bool_t result = kTRUE;
  Float_t kEpsilon=0.00000000000000001;
  
  // calcul the mean value per col
  for(Int_t icol = 0; icol < fNcols; icol++){

    Float_t mean = 0.0;
    Float_t nb   = 0.0;

    for(Int_t irow = 0; irow < fNrows; irow++){
      if((GetValue(icol,irow) > 0.06) && (GetValue(icol,irow) < 0.15)){
	mean += GetValue(icol,irow);
	nb += 1.0;
      }
    }
    
    if(nb > kEpsilon) {
      
      mean = mean/nb;

      if(mean > kEpsilon){
	for(Int_t irow = 0; irow < fNrows; irow++){
	  Float_t value = GetValue(icol,irow);
	  SetValue(icol,irow,(Float_t)(value/mean));
	}    
      }
      else result = kFALSE;
      
    }

    else result = kFALSE;

    
  }
 

  return result;
}

//__________________________________________________________________________________
TH2F * AliTRDCalROC::MakeHisto2D(Float_t min, Float_t max,Int_t type,  Float_t mu)
{
  //
  // make 2D histo
  // type -1 = user defined range
  //       0 = nsigma cut nsigma=min
  //       1 = delta cut around median delta=min

  Float_t kEpsilonr = 0.005;
  gStyle->SetPalette(1);
  
  if (type>=0){
    if (type==0){
      // nsigma range
      Float_t mean  = GetMean();
      Float_t sigma = GetRMS();
      Float_t nsigma = TMath::Abs(min);
      sigma *= nsigma;
      if(sigma < kEpsilonr) sigma = kEpsilonr;
      min = mean-sigma;
      max = mean+sigma;
    }
    if (type==1){
      // fixed range
      Float_t mean   = GetMedian();
      Float_t  delta = min;
      min = mean-delta;
      max = mean+delta;
    }
    if (type==2){
      Double_t sigma;
      Float_t mean  = GetLTM(&sigma,max);
      sigma*=min;
      if(sigma < kEpsilonr) sigma = kEpsilonr;
      min = mean-sigma;
      max = mean+sigma;

    }
  }
  
  char  name[1000];
  sprintf(name,"%s 2D Plane %d Chamber %d",GetTitle(),fPla, fCha);
  TH2F * his = new TH2F(name,name,fNrows,0, fNrows, fNcols, 0, fNcols);
  for (Int_t irow=0; irow<fNrows; irow++){
    for (Int_t icol=0; icol<fNcols; icol++){
      his->Fill(irow+0.5,icol+0.5,GetValue(icol,irow)*mu);
    }
  }
  his->SetStats(0);
  his->SetMaximum(max);
  his->SetMinimum(min);
  return his;
}

//_______________________________________________________________________________________
TH1F * AliTRDCalROC::MakeHisto1D(Float_t min, Float_t max,Int_t type,  Float_t mu)
{
  //
  // make 1D histo
  // type -1 = user defined range
  //       0 = nsigma cut nsigma=min
  //       1 = delta cut around median delta=min

  Float_t kEpsilonr = 0.005;

  if (type>=0){
    if (type==0){
      // nsigma range
      Float_t mean  = GetMean();
      Float_t sigma = GetRMS();
      Float_t nsigma = TMath::Abs(min);
      sigma *= nsigma;
      if(sigma < kEpsilonr) sigma = kEpsilonr;
      min = mean-sigma;
      max = mean+sigma;
    }
    if (type==1){
      // fixed range
      Float_t mean   = GetMedian();
      Float_t  delta = min;
      min = mean-delta;
      max = mean+delta;
    }
    if (type==2){
      //
      // LTM mean +- nsigma
      //
      Double_t sigma;
      Float_t mean  = GetLTM(&sigma,max);
      sigma*=min;
      if(sigma < kEpsilonr) sigma = kEpsilonr;
      min = mean-sigma;
      max = mean+sigma;
    }
  }
  char  name[1000];
  sprintf(name,"%s 1D Plane %d Chamber %d",GetTitle(),fPla, fCha);
  TH1F * his = new TH1F(name,name,100, min,max);
  for (Int_t irow=0; irow<fNrows; irow++){
    for (Int_t icol=0; icol<fNcols; icol++){
      his->Fill(GetValue(icol,irow)*mu);
    }
  }
  return his;
}
Commit	Line	Data
7754cd1f	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	///////////////////////////////////////////////////////////////////////////////
	19	// //
	20	// Calibration base class for a single ROC //
	21	// Contains one UShort_t value per pad //
	22	// However, values are set and get as float, there are stored internally as //
	23	// (UShort_t) value * 10000 //
	24	// //
	25	///////////////////////////////////////////////////////////////////////////////
	26
4c87f26e	27	#include <TStyle.h>
c893147d	28	#include <TMath.h>
	29	#include <TH2F.h>
	30	#include <TH1F.h>
4c87f26e	31
4c87f26e	32	#include "AliMathBase.h"
c893147d	33
c893147d	34	#include "AliTRDCalROC.h"
7754cd1f	35
	36	ClassImp(AliTRDCalROC)
	37
	38	//_____________________________________________________________________________
2745a409	39	AliTRDCalROC::AliTRDCalROC()
	40	:TObject()
	41	,fPla(0)
	42	,fCha(0)
	43	,fNrows(0)
	44	,fNcols(0)
	45	,fNchannels(0)
	46	,fData(0)
7754cd1f	47	{
	48	//
	49	// Default constructor
	50	//
	51
7754cd1f	52	}
	53
	54	//_____________________________________________________________________________
2745a409	55	AliTRDCalROC::AliTRDCalROC(Int_t p, Int_t c)
	56	:TObject()
	57	,fPla(p)
	58	,fCha(c)
	59	,fNrows(0)
	60	,fNcols(144)
	61	,fNchannels(0)
	62	,fData(0)
7754cd1f	63	{
	64	//
	65	// Constructor that initializes a given pad plane type
	66	//
	67
7754cd1f	68	//
	69	// The pad plane parameter
	70	//
	71	switch (p) {
	72	case 0:
	73	if (c == 2) {
	74	// L0C0 type
	75	fNrows = 12;
	76	}
	77	else {
	78	// L0C1 type
	79	fNrows = 16;
	80	}
	81	break;
	82	case 1:
	83	if (c == 2) {
	84	// L1C0 type
	85	fNrows = 12;
	86	}
	87	else {
	88	// L1C1 type
	89	fNrows = 16;
	90	}
	91	break;
	92	case 2:
	93	if (c == 2) {
	94	// L2C0 type
	95	fNrows = 12;
	96	}
	97	else {
	98	// L2C1 type
	99	fNrows = 16;
	100	}
	101	break;
	102	case 3:
	103	if (c == 2) {
	104	// L3C0 type
	105	fNrows = 12;
	106	}
	107	else {
	108	// L3C1 type
	109	fNrows = 16;
	110	}
	111	break;
	112	case 4:
	113	if (c == 2) {
	114	// L4C0 type
	115	fNrows = 12;
	116	}
	117	else {
	118	// L4C1 type
	119	fNrows = 16;
	120	}
	121	break;
	122	case 5:
	123	if (c == 2) {
	124	// L5C0 type
	125	fNrows = 12;
	126	}
	127	else {
	128	// L5C1 type
	129	fNrows = 16;
	130	}
	131	break;
132	};
133
134	fNchannels = fNrows * fNcols;
2745a409	135	if (fNchannels != 0) {
7754cd1f	136	fData = new UShort_t[fNchannels];
2745a409	137	}
7754cd1f	138
2745a409	139	for (Int_t i=0; i<fNchannels; ++i) {
7754cd1f	140	fData[i] = 0;
2745a409	141	}
2745a409	142
7754cd1f	143	}
	144
	145	//_____________________________________________________________________________
2745a409	146	AliTRDCalROC::AliTRDCalROC(const AliTRDCalROC &c)
	147	:TObject(c)
	148	,fPla(c.fPla)
	149	,fCha(c.fCha)
	150	,fNrows(c.fNrows)
	151	,fNcols(c.fNcols)
	152	,fNchannels(c.fNchannels)
	153	,fData(0)
7754cd1f	154	{
	155	//
	156	// AliTRDCalROC copy constructor
	157	//
	158
2745a409	159	Int_t iBin = 0;
2745a409	160
4c87f26e	161	fData = new UShort_t[fNchannels];
2745a409	162	for (iBin = 0; iBin < fNchannels; iBin++) {
4c87f26e	163	fData[iBin] = ((AliTRDCalROC &) c).fData[iBin];
2745a409	164	}
7754cd1f	165
	166	}
	167
	168	//_____________________________________________________________________________
	169	AliTRDCalROC::~AliTRDCalROC()
	170	{
	171	//
	172	// AliTRDCalROC destructor
	173	//
	174
	175	if (fData) {
	176	delete [] fData;
	177	fData = 0;
	178	}
2745a409	179
7754cd1f	180	}
	181
	182	//_____________________________________________________________________________
	183	AliTRDCalROC &AliTRDCalROC::operator=(const AliTRDCalROC &c)
	184	{
	185	//
	186	// Assignment operator
	187	//
	188
	189	if (this != &c) ((AliTRDCalROC &) c).Copy(*this);
	190	return *this;
	191
	192	}
	193
	194	//_____________________________________________________________________________
	195	void AliTRDCalROC::Copy(TObject &c) const
	196	{
	197	//
	198	// Copy function
	199	//
	200
	201	((AliTRDCalROC &) c).fPla = fPla;
	202	((AliTRDCalROC &) c).fCha = fCha;
	203
	204	((AliTRDCalROC &) c).fNrows = fNrows;
	205	((AliTRDCalROC &) c).fNcols = fNcols;
	206
	207	Int_t iBin = 0;
	208
	209	((AliTRDCalROC &) c).fNchannels = fNchannels;
	210
	211	if (((AliTRDCalROC &) c).fData) delete [] ((AliTRDCalROC &) c).fData;
	212	((AliTRDCalROC &) c).fData = new UShort_t[fNchannels];
	213	for (iBin = 0; iBin < fNchannels; iBin++) {
	214	((AliTRDCalROC &) c).fData[iBin] = fData[iBin];
	215	}
4c87f26e	216
7754cd1f	217	TObject::Copy(c);
	218
	219	}
	220
4c87f26e	221	//___________________________________________________________________________________
c893147d	222	Double_t AliTRDCalROC::GetMean(AliTRDCalROC* const outlierROC) const
4c87f26e	223	{
	224	//
	225	// Calculate the mean
	226	//
	227
	228	Double_t *ddata = new Double_t[fNchannels];
c893147d	229	Int_t nPoints = 0;
4c87f26e	230	for (Int_t i=0;i<fNchannels;i++) {
4c87f26e	231	if ((!outlierROC) \|\| (!(outlierROC->GetValue(i)))) {
37b0cf5e	232	if(fData[i] > 0.000000000000001){
c893147d	233	ddata[nPoints]= (Double_t) fData[i]/10000;
c893147d	234	nPoints++;
37b0cf5e	235	}
4c87f26e	236	}
4c87f26e	237	}
c893147d	238	Double_t mean = TMath::Mean(nPoints,ddata);
4c87f26e	239	delete [] ddata;
	240	return mean;
	241	}
	242
	243	//_______________________________________________________________________________________
c893147d	244	Double_t AliTRDCalROC::GetMedian(AliTRDCalROC* const outlierROC) const
4c87f26e	245	{
	246	//
	247	// Calculate the median
	248	//
	249
	250	Double_t *ddata = new Double_t[fNchannels];
c893147d	251	Int_t nPoints = 0;
4c87f26e	252	for (Int_t i=0;i<fNchannels;i++) {
4c87f26e	253	if ((!outlierROC) \|\| (!(outlierROC->GetValue(i)))) {
37b0cf5e	254	if(fData[i] > 0.000000000000001){
c893147d	255	ddata[nPoints]= (Double_t) fData[i]/10000;
c893147d	256	nPoints++;
37b0cf5e	257	}
37b0cf5e	258	}
4c87f26e	259	}
c893147d	260	Double_t mean = TMath::Median(nPoints,ddata);
4c87f26e	261	delete [] ddata;
	262	return mean;
	263	}
	264
	265	//____________________________________________________________________________________________
c893147d	266	Double_t AliTRDCalROC::GetRMS(AliTRDCalROC* const outlierROC) const
4c87f26e	267	{
	268	//
	269	// Calculate the RMS
	270	//
	271
	272	Double_t *ddata = new Double_t[fNchannels];
c893147d	273	Int_t nPoints = 0;
37b0cf5e	274	for (Int_t i=0;i<fNchannels;i++) {
	275	if ((!outlierROC) \|\| (!(outlierROC->GetValue(i)))) {
	276	if(fData[i] > 0.000000000000001){
c893147d	277	ddata[nPoints]= (Double_t)fData[i]/10000;
c893147d	278	nPoints++;
37b0cf5e	279	}
	280	}
	281	}
c893147d	282	Double_t mean = TMath::RMS(nPoints,ddata);
37b0cf5e	283	delete [] ddata;
37b0cf5e	284	return mean;
4c87f26e	285	}
	286
	287	//______________________________________________________________________________________________
c893147d	288	Double_t AliTRDCalROC::GetLTM(Double_t sigma, Double_t fraction, AliTRDCalROC const outlierROC)
7754cd1f	289	{
7754cd1f	290	//
4c87f26e	291	// Calculate LTM mean and sigma
7754cd1f	292	//
7754cd1f	293
4c87f26e	294	Double_t *ddata = new Double_t[fNchannels];
4c87f26e	295	Double_t mean=0, lsigma=0;
c893147d	296	UInt_t nPoints = 0;
4c87f26e	297	for (Int_t i=0;i<fNchannels;i++) {
4c87f26e	298	if (!outlierROC \|\| !(outlierROC->GetValue(i))) {
37b0cf5e	299	if(fData[i] > 0.000000000000001){
c893147d	300	ddata[nPoints]= (Double_t) fData[i]/10000;
c893147d	301	nPoints++;
37b0cf5e	302	}
4c87f26e	303	}
7754cd1f	304	}
c893147d	305	Int_t hh = TMath::Min(TMath::Nint(fraction *nPoints), Int_t(nPoints));
c893147d	306	AliMathBase::EvaluateUni(nPoints,ddata, mean, lsigma, hh);
4c87f26e	307	if (sigma) *sigma=lsigma;
	308	delete [] ddata;
	309	return mean;
	310	}
2745a409	311
4c87f26e	312	//___________________________________________________________________________________
	313	Bool_t AliTRDCalROC::Add(Float_t c1)
	314	{
	315	//
	316	// add constant
	317	//
c893147d	318
4c87f26e	319	Bool_t result = kTRUE;
	320	for (Int_t idata = 0; idata< fNchannels; idata++) {
	321	if(((GetValue(idata)+c1) <= 6.5535) && ((GetValue(idata)+c1) >= 0.0)) SetValue(idata,GetValue(idata)+c1);
	322	else {
	323	SetValue(idata,0.0);
	324	result = kFALSE;
	325	}
	326	}
	327	return result;
	328	}
	329
	330	//_______________________________________________________________________________________
	331	Bool_t AliTRDCalROC::Multiply(Float_t c1)
	332	{
	333	//
	334	// multiply constant
	335	//
c893147d	336
4c87f26e	337	Bool_t result = kTRUE;
	338	if(c1 < 0) return kFALSE;
	339	for (Int_t idata = 0; idata< fNchannels; idata++) {
	340	if((GetValue(idata)c1) <= 6.5535) SetValue(idata,GetValue(idata)c1);
	341	else {
	342	SetValue(idata,0.0);
	343	result = kFALSE;
	344	}
	345	}
	346	return result;
	347	}
	348
	349	//____________________________________________________________________________________________
	350	Bool_t AliTRDCalROC::Add(const AliTRDCalROC * roc, Double_t c1)
	351	{
	352	//
	353	// add values
	354	//
c893147d	355
4c87f26e	356	Bool_t result = kTRUE;
	357	for (Int_t idata = 0; idata< fNchannels; idata++){
	358	if(((GetValue(idata)+roc->GetValue(idata)*c1) <= 6.5535) &&
	359	((GetValue(idata)+roc->GetValue(idata)*c1) >= 0.0))
	360	SetValue(idata,GetValue(idata)+roc->GetValue(idata)*c1);
	361	else {
	362	SetValue(idata,0.0);
	363	result = kFALSE;
	364	}
	365	}
	366	return result;
	367	}
	368
	369	//____________________________________________________________________________________________
	370	Bool_t AliTRDCalROC::Multiply(const AliTRDCalROC* roc)
	371	{
	372	//
	373	// multiply values - per by pad
	374	//
c893147d	375
4c87f26e	376	Bool_t result = kTRUE;
	377	for (Int_t idata = 0; idata< fNchannels; idata++){
	378	if((GetValue(idata)*roc->GetValue(idata)) <= 6.5535)
	379	SetValue(idata,GetValue(idata)*roc->GetValue(idata));
	380	else {
	381	SetValue(idata,0.0);
	382	result = kFALSE;
	383	}
	384	}
	385	return result;
	386	}
	387
	388	//______________________________________________________________________________________________
	389	Bool_t AliTRDCalROC::Divide(const AliTRDCalROC* roc)
	390	{
	391	//
	392	// divide values
	393	//
c893147d	394
4c87f26e	395	Bool_t result = kTRUE;
	396	Float_t kEpsilon=0.00000000000000001;
	397	for (Int_t idata = 0; idata< fNchannels; idata++){
	398	if (TMath::Abs(roc->GetValue(idata))>kEpsilon){
	399	if((GetValue(idata)/roc->GetValue(idata)) <= 6.5535)
	400	SetValue(idata,GetValue(idata)/roc->GetValue(idata));
	401	else {
	402	SetValue(idata,0.0);
	403	result = kFALSE;
	404	}
	405	}
	406	else result = kFALSE;
	407	}
	408	return result;
	409	}
37b0cf5e	410	//______________________________________________________________________________________________
	411	Bool_t AliTRDCalROC::Unfold()
	412	{
	413	//
	414	// Compute the mean value per pad col
	415	// Divide with this value each pad col
	416	// This is for the noise study
	417	// Return kFALSE if one or more of the pad col was not normalised
	418	//
c893147d	419
37b0cf5e	420	Bool_t result = kTRUE;
	421	Float_t kEpsilon=0.00000000000000001;
	422
	423	// calcul the mean value per col
	424	for(Int_t icol = 0; icol < fNcols; icol++){
	425
	426	Float_t mean = 0.0;
	427	Float_t nb = 0.0;
	428
	429	for(Int_t irow = 0; irow < fNrows; irow++){
	430	if((GetValue(icol,irow) > 0.06) && (GetValue(icol,irow) < 0.15)){
	431	mean += GetValue(icol,irow);
	432	nb += 1.0;
	433	}
	434	}
	435
	436	if(nb > kEpsilon) {
	437
	438	mean = mean/nb;
	439
	440	if(mean > kEpsilon){
	441	for(Int_t irow = 0; irow < fNrows; irow++){
	442	Float_t value = GetValue(icol,irow);
	443	SetValue(icol,irow,(Float_t)(value/mean));
	444	}
	445	}
	446	else result = kFALSE;
	447
	448	}
	449
	450	else result = kFALSE;
	451
	452
	453	}
	454
	455
	456	return result;
	457	}
4c87f26e	458
4c87f26e	459	//__________________________________________________________________________________
7dcf6933	460	TH2F * AliTRDCalROC::MakeHisto2D(Float_t min, Float_t max,Int_t type, Float_t mu)
4c87f26e	461	{
	462	//
	463	// make 2D histo
	464	// type -1 = user defined range
	465	// 0 = nsigma cut nsigma=min
	466	// 1 = delta cut around median delta=min
c893147d	467
4c87f26e	468	Float_t kEpsilonr = 0.005;
	469	gStyle->SetPalette(1);
	470
	471	if (type>=0){
	472	if (type==0){
	473	// nsigma range
	474	Float_t mean = GetMean();
	475	Float_t sigma = GetRMS();
	476	Float_t nsigma = TMath::Abs(min);
	477	sigma *= nsigma;
	478	if(sigma < kEpsilonr) sigma = kEpsilonr;
	479	min = mean-sigma;
	480	max = mean+sigma;
	481	}
	482	if (type==1){
	483	// fixed range
	484	Float_t mean = GetMedian();
	485	Float_t delta = min;
	486	min = mean-delta;
	487	max = mean+delta;
	488	}
	489	if (type==2){
	490	Double_t sigma;
	491	Float_t mean = GetLTM(&sigma,max);
	492	sigma*=min;
	493	if(sigma < kEpsilonr) sigma = kEpsilonr;
	494	min = mean-sigma;
	495	max = mean+sigma;
	496
	497	}
	498	}
	499
	500	char name[1000];
	501	sprintf(name,"%s 2D Plane %d Chamber %d",GetTitle(),fPla, fCha);
	502	TH2F * his = new TH2F(name,name,fNrows,0, fNrows, fNcols, 0, fNcols);
	503	for (Int_t irow=0; irow<fNrows; irow++){
	504	for (Int_t icol=0; icol<fNcols; icol++){
7dcf6933	505	his->Fill(irow+0.5,icol+0.5,GetValue(icol,irow)*mu);
4c87f26e	506	}
	507	}
	508	his->SetStats(0);
	509	his->SetMaximum(max);
	510	his->SetMinimum(min);
	511	return his;
	512	}
	513
	514	//_______________________________________________________________________________________
7dcf6933	515	TH1F * AliTRDCalROC::MakeHisto1D(Float_t min, Float_t max,Int_t type, Float_t mu)
4c87f26e	516	{
	517	//
	518	// make 1D histo
	519	// type -1 = user defined range
	520	// 0 = nsigma cut nsigma=min
	521	// 1 = delta cut around median delta=min
4c87f26e	522
c893147d	523	Float_t kEpsilonr = 0.005;
4c87f26e	524
	525	if (type>=0){
	526	if (type==0){
	527	// nsigma range
	528	Float_t mean = GetMean();
	529	Float_t sigma = GetRMS();
	530	Float_t nsigma = TMath::Abs(min);
	531	sigma *= nsigma;
	532	if(sigma < kEpsilonr) sigma = kEpsilonr;
	533	min = mean-sigma;
	534	max = mean+sigma;
	535	}
	536	if (type==1){
	537	// fixed range
	538	Float_t mean = GetMedian();
	539	Float_t delta = min;
	540	min = mean-delta;
	541	max = mean+delta;
	542	}
	543	if (type==2){
	544	//
	545	// LTM mean +- nsigma
	546	//
	547	Double_t sigma;
	548	Float_t mean = GetLTM(&sigma,max);
	549	sigma*=min;
	550	if(sigma < kEpsilonr) sigma = kEpsilonr;
	551	min = mean-sigma;
	552	max = mean+sigma;
	553	}
	554	}
	555	char name[1000];
	556	sprintf(name,"%s 1D Plane %d Chamber %d",GetTitle(),fPla, fCha);
	557	TH1F * his = new TH1F(name,name,100, min,max);
	558	for (Int_t irow=0; irow<fNrows; irow++){
	559	for (Int_t icol=0; icol<fNcols; icol++){
7dcf6933	560	his->Fill(GetValue(icol,irow)*mu);
4c87f26e	561	}
	562	}
	563	return his;
7754cd1f	564	}