// //
///////////////////////////////////////////////////////////////////////////////
+#include <TMath.h>
+#include <TH1F.h>
+#include <TH2F.h>
+#include <TStyle.h>
+
#include "AliTRDCalDet.h"
+#include "AliTRDgeometry.h"
+#include "AliMathBase.h"
+#include "AliTRDpadPlane.h"
ClassImp(AliTRDCalDet)
}
+//___________________________________________________________________________________
+Double_t AliTRDCalDet::GetMean(AliTRDCalDet* outlierDet)
+{
+ //
+ // Calculate the mean
+ //
+
+ if (!outlierDet) return TMath::Mean(kNdet,fData);
+ Double_t *ddata = new Double_t[kNdet];
+ Int_t NPoints = 0;
+ for (Int_t i=0;i<kNdet;i++) {
+ if (!(outlierDet->GetValue(i))) {
+ ddata[NPoints]= fData[NPoints];
+ NPoints++;
+ }
+ }
+ Double_t mean = TMath::Mean(NPoints,ddata);
+ delete [] ddata;
+ return mean;
+}
+
+//_______________________________________________________________________________________
+Double_t AliTRDCalDet::GetMedian(AliTRDCalDet* outlierDet)
+{
+ //
+ // Calculate the median
+ //
+
+ if (!outlierDet) return (Double_t) TMath::Median(kNdet,fData);
+ Double_t *ddata = new Double_t[kNdet];
+ Int_t NPoints = 0;
+ for (Int_t i=0;i<kNdet;i++) {
+ if (!(outlierDet->GetValue(i))) {
+ ddata[NPoints]= fData[NPoints];
+ NPoints++;
+ }
+ }
+ Double_t mean = TMath::Median(NPoints,ddata);
+ delete [] ddata;
+ return mean;
+
+}
+
+//____________________________________________________________________________________________
+Double_t AliTRDCalDet::GetRMS(AliTRDCalDet* outlierDet)
+{
+ //
+ // Calculate the RMS
+ //
+
+ if (!outlierDet) return TMath::RMS(kNdet,fData);
+ Double_t *ddata = new Double_t[kNdet];
+ Int_t NPoints = 0;
+ for (Int_t i=0;i<kNdet;i++) {
+ if (!(outlierDet->GetValue(i))) {
+ ddata[NPoints]= fData[NPoints];
+ NPoints++;
+ }
+ }
+ Double_t mean = TMath::RMS(NPoints,ddata);
+ delete [] ddata;
+ return mean;
+}
+
+//______________________________________________________________________________________________
+Double_t AliTRDCalDet::GetLTM(Double_t *sigma, Double_t fraction, AliTRDCalDet* outlierDet)
+{
+ //
+ // Calculate LTM mean and sigma
+ //
+
+ Double_t *ddata = new Double_t[kNdet];
+ Double_t mean=0, lsigma=0;
+ UInt_t NPoints = 0;
+ for (Int_t i=0;i<kNdet;i++) {
+ if (!outlierDet || !(outlierDet->GetValue(i))) {
+ ddata[NPoints]= fData[NPoints];
+ 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;
+}
+
+//_________________________________________________________________________
+TH1F * AliTRDCalDet::MakeHisto1Distribution(Float_t min, Float_t max,Int_t type)
+{
+ //
+ // make 1D histo
+ // type -1 = user defined range
+ // 0 = nsigma cut nsigma=min
+ // 1 = delta cut around median delta=min
+ //
+
+ if (type>=0){
+ if (type==0){
+ // nsigma range
+ Float_t mean = GetMean();
+ Float_t sigma = GetRMS();
+ Float_t nsigma = TMath::Abs(min);
+ min = mean-nsigma*sigma;
+ max = mean+nsigma*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;
+ min = mean-sigma;
+ max = mean+sigma;
+ }
+ }
+ char name[1000];
+ sprintf(name,"%s CalDet 1Distribution",GetTitle());
+ TH1F * his = new TH1F(name,name,100, min,max);
+ for (Int_t idet=0; idet<kNdet; idet++){
+ his->Fill(GetValue(idet));
+ }
+ return his;
+}
+
+//________________________________________________________________________________
+TH1F * AliTRDCalDet::MakeHisto1DAsFunctionOfDet(Float_t min, Float_t max,Int_t type)
+{
+ //
+ // make 1D histo
+ // type -1 = user defined range
+ // 0 = nsigma cut nsigma=min
+ // 1 = delta cut around median delta=min
+ //
+
+ if (type>=0){
+ if (type==0){
+ // nsigma range
+ Float_t mean = GetMean();
+ Float_t sigma = GetRMS();
+ Float_t nsigma = TMath::Abs(min);
+ min = mean-nsigma*sigma;
+ max = mean+nsigma*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;
+ min = mean-sigma;
+ max = mean+sigma;
+
+ }
+ }
+
+ char name[1000];
+ sprintf(name,"%s CalDet as function of det",GetTitle());
+ TH1F * his = new TH1F(name,name,kNdet, 0, kNdet);
+ for(Int_t det = 0; det< kNdet; det++){
+ his->Fill(det+0.5,GetValue(det));
+ }
+
+ his->SetMaximum(max);
+ his->SetMinimum(min);
+ return his;
+}
+
+//_____________________________________________________________________________
+TH2F *AliTRDCalDet::MakeHisto2DCh(Int_t ch, Float_t min, Float_t max, Int_t type)
+{
+ //
+ // Make 2D graph
+ // ch - chamber
+ // type -1 = user defined range
+ // 0 = nsigma cut nsigma=min
+ // 1 = delta cut around median delta=min
+ //
+
+ 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);
+ min = mean-nsigma*sigma;
+ max = mean+nsigma*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;
+ min = mean-sigma;
+ max = mean+sigma;
+
+ }
+ }
+
+ AliTRDgeometry *TRDgeo = new AliTRDgeometry();
+ TRDgeo->Init();
+
+
+ Double_t poslocal[3] = {0.0,0.0,0.0};
+ Double_t posglobal[3] = {0.0,0.0,0.0};
+
+ char name[1000];
+ sprintf(name,"%s CalDet 2D ch %d",GetTitle(),ch);
+ TH2F * his = new TH2F(name, name, 400,-400.0,400.0,400,-400.0,400.0);
+
+ // Where we begin
+ Int_t offsetch = 6*ch;
+
+
+ for (Int_t isec = 0; isec < kNsect; isec++){
+ for(Int_t ipl = 0; ipl < kNplan; ipl++){
+ Int_t det = offsetch+isec*30+ipl;
+ AliTRDpadPlane *padPlane = new AliTRDpadPlane(ipl,ch);
+ for (Int_t icol=0; icol<padPlane->GetNcols(); icol++){
+ poslocal[0] = TRDgeo->GetTime0(ipl);
+ poslocal[2] = padPlane->GetRowPos(0);
+ poslocal[1] = padPlane->GetColPos(icol);
+ TRDgeo->RotateBack(det,poslocal,posglobal);
+ Int_t binx = 1+TMath::Nint((posglobal[0]+400.0)*0.5);
+ Int_t biny = 1+TMath::Nint((posglobal[1]+400.0)*0.5);
+ his->SetBinContent(binx,biny,fData[det]);
+ }
+ }
+ }
+ his->SetXTitle("x (cm)");
+ his->SetYTitle("y (cm)");
+ his->SetStats(0);
+ his->SetMaximum(max);
+ his->SetMinimum(min);
+ return his;
+}
+
+//_____________________________________________________________________________
+TH2F *AliTRDCalDet::MakeHisto2DSmPl(Int_t sm, Int_t pl, Float_t min, Float_t max, Int_t type)
+{
+ //
+ // Make 2D graph
+ // sm - supermodule number
+ // pl - plane number
+ // type -1 = user defined range
+ // 0 = nsigma cut nsigma=min
+ // 1 = delta cut around median delta=min
+ //
+
+ 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);
+ min = mean-nsigma*sigma;
+ max = mean+nsigma*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;
+ min = mean-sigma;
+ max = mean+sigma;
+
+ }
+ }
+
+ AliTRDpadPlane *padPlane0 = new AliTRDpadPlane(pl,0);
+ Double_t row0 = padPlane0->GetRow0();
+ Double_t col0 = padPlane0->GetCol0();
+
+
+ char name[1000];
+ sprintf(name,"%s CalDet 2D sm %d and pl %d",GetTitle(),sm,pl);
+ TH2F * his = new TH2F( name, name, 5, -TMath::Abs(row0), TMath::Abs(row0)
+ , 4,-2*TMath::Abs(col0),2*TMath::Abs(col0));
+
+ // Where we begin
+ Int_t offsetsmpl = 30*sm+pl;
+
+
+ for (Int_t k = 0; k < kNcham; k++){
+ Int_t det = offsetsmpl+k*6;
+ Int_t kb = kNcham-1-k;
+ his->SetBinContent(kb+1,2,fData[det]);
+ his->SetBinContent(kb+1,3,fData[det]);
+ }
+ his->SetXTitle("z (cm)");
+ his->SetYTitle("y (cm)");
+ his->SetStats(0);
+ his->SetMaximum(max);
+ his->SetMinimum(min);
+ return his;
+}
+
+//_____________________________________________________________________________
+void AliTRDCalDet::Add(Float_t c1)
+{
+ //
+ // Add constant for all detectors
+ //
+
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ fData[idet] += c1;
+ }
+}
+
+//_____________________________________________________________________________
+void AliTRDCalDet::Multiply(Float_t c1)
+{
+ //
+ // multiply constant for all detectors
+ //
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ fData[idet] *= c1;
+ }
+}
+
+//_____________________________________________________________________________
+void AliTRDCalDet::Add(const AliTRDCalDet * calDet, Double_t c1)
+{
+ //
+ // add caldet channel by channel multiplied by c1
+ //
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ fData[idet] += calDet->GetValue(idet)*c1;
+ }
+}
+
+//_____________________________________________________________________________
+void AliTRDCalDet::Multiply(const AliTRDCalDet * calDet)
+{
+ //
+ // multiply caldet channel by channel
+ //
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ fData[idet] *= calDet->GetValue(idet);
+ }
+}
+
+//_____________________________________________________________________________
+void AliTRDCalDet::Divide(const AliTRDCalDet * calDet)
+{
+ //
+ // divide caldet channel by channel
+ //
+ Float_t kEpsilon=0.00000000000000001;
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if (TMath::Abs(calDet->GetValue(idet))>kEpsilon){
+ fData[idet] /= calDet->GetValue(idet);
+ }
+ }
+}
///////////////////////////////////////////////////////////////////////////////
#include "TNamed.h"
-#include "AliTRDgeometry.h"
+#include "TMath.h"
+#include "../AliTRDgeometry.h"
+
+class TH1F;
+class TH2F;
class AliTRDCalDet : public TNamed {
void SetValue(Int_t d, Float_t value) { fData[d] = value; };
void SetValue(Int_t p, Int_t c, Int_t s, Float_t value)
{ fData[AliTRDgeometry::GetDetector(p,c,s)] = value; };
-
+
+ // statistic
+ Double_t GetMean(AliTRDCalDet *outlierDet=0);
+ Double_t GetRMS(AliTRDCalDet *outlierDet=0);
+ Double_t GetMedian(AliTRDCalDet *outlierDet=0);
+ Double_t GetLTM(Double_t *sigma=0, Double_t fraction=0.9, AliTRDCalDet *outlierDet=0);
+
+ // Plot functions
+ TH1F * MakeHisto1Distribution(Float_t min=4, Float_t max=-4, Int_t type=0);
+ TH1F * MakeHisto1DAsFunctionOfDet(Float_t min=4, Float_t max=-4, Int_t type=0);
+ TH2F * MakeHisto2DCh(Int_t ch, Float_t min=4, Float_t max=-4, Int_t type=0);
+ TH2F * MakeHisto2DSmPl(Int_t sm, Int_t pl, Float_t min=4, Float_t max=-4, Int_t type=0);
+
+ // algebra functions
+ void Add(Float_t c1);
+ void Multiply(Float_t c1);
+ void Add(const AliTRDCalDet * calDet, Double_t c1 = 1);
+ void Multiply(const AliTRDCalDet * calDet);
+ void Divide(const AliTRDCalDet * calDet);
+
protected:
Float_t fData[kNdet]; //[kNdet] Data
// //
///////////////////////////////////////////////////////////////////////////////
+#include <TMath.h>
+#include <TH2F.h>
+#include <TH1F.h>
+#include <TStyle.h>
+
#include "AliTRDCalPad.h"
#include "AliTRDCalROC.h"
#include "AliTRDCalDet.h"
+#include "AliTRDpadPlane.h"
+#include "AliMathBase.h"
+#include "AliTRDgeometry.h"
ClassImp(AliTRDCalPad)
}
//_____________________________________________________________________________
-AliTRDCalPad::AliTRDCalPad(const AliTRDCalPad &c):TNamed(c)
+AliTRDCalPad::AliTRDCalPad(const AliTRDCalPad &c)
+ :TNamed(c)
{
//
// AliTRDCalPad copy constructor
//
- ((AliTRDCalPad &) c).Copy(*this);
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ fROC[idet] = new AliTRDCalROC(*((AliTRDCalPad &) c).fROC[idet]);
+ }
}
-///_____________________________________________________________________________
+//_____________________________________________________________________________
AliTRDCalPad::~AliTRDCalPad()
{
//
}
//_____________________________________________________________________________
-void AliTRDCalPad::ScaleROCs(AliTRDCalDet* values)
+Bool_t AliTRDCalPad::ScaleROCs(const AliTRDCalDet* values)
{
//
// Scales ROCs of this class with the values from the class <values>
//
if (!values)
- return;
-
+ return kFALSE;
+ Bool_t result = kTRUE;
for (Int_t idet = 0; idet < kNdet; idet++) {
if (fROC[idet]) {
- fROC[idet]->Scale(values->GetValue(idet));
+ if(!fROC[idet]->Multiply(values->GetValue(idet))) result = kFALSE;
+ }
+ }
+ return result;
+}
+
+//_____________________________________________________________________________
+Double_t AliTRDCalPad::GetMeanRMS(Double_t &rms, const AliTRDCalDet *calDet, Int_t type)
+{
+ //
+ // Calculate mean an RMS of all rocs
+ // If calDet correct the CalROC from the detector coefficient
+ // type == 0 for gain and vdrift
+ // type == 1 for t0
+ //
+ Double_t factor = 0.0;
+ if(type == 0) factor = 1.0;
+ Double_t sum = 0, sum2 = 0, n=0, val=0;
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet) factor = calDet->GetValue(idet);
+ AliTRDCalROC *calRoc = fROC[idet];
+ if ( calRoc ){
+ for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
+ for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
+ if(type == 0) val = calRoc->GetValue(icol,irow)*factor;
+ else val = calRoc->GetValue(icol,irow)+factor;
+ sum+=val;
+ sum2+=val*val;
+ n++;
+ }
+ }
+ }
+ }
+ Double_t n1 = 1./n;
+ Double_t mean = sum*n1;
+ rms = TMath::Sqrt(TMath::Abs(sum2*n1-mean*mean));
+ return mean;
+}
+
+//_____________________________________________________________________________
+Double_t AliTRDCalPad::GetMean(const AliTRDCalDet *calDet, Int_t type, AliTRDCalPad* outlierPad)
+{
+ //
+ // return mean of the mean of all ROCs
+ // If calDet correct the CalROC from the detector coefficient
+ // type == 0 for gain and vdrift
+ // type == 1 for t0
+ //
+ Double_t factor = 0.0;
+ if(type == 0) factor = 1.0;
+ Double_t arr[kNdet];
+ Int_t n=0;
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet) factor = calDet->GetValue(idet);
+ AliTRDCalROC *calRoc = fROC[idet];
+ if ( calRoc ){
+ AliTRDCalROC* outlierROC = 0;
+ if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
+ if(type == 0) arr[n] = calRoc->GetMean(outlierROC)*factor;
+ else arr[n] = calRoc->GetMean(outlierROC)+factor;
+ n++;
+ }
+ }
+ return TMath::Mean(n,arr);
+}
+
+//_____________________________________________________________________________
+Double_t AliTRDCalPad::GetRMS(const AliTRDCalDet *calDet, Int_t type, AliTRDCalPad* outlierPad)
+{
+ //
+ // return mean of the RMS of all ROCs
+ // If calDet correct the CalROC from the detector coefficient
+ // type == 0 for gain and vdrift
+ // type == 1 for t0
+ //
+ Double_t factor = 0.0;
+ if(type == 0) factor = 1.0;
+ Double_t arr[kNdet];
+ Int_t n=0;
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet) factor = calDet->GetValue(idet);
+ AliTRDCalROC *calRoc = fROC[idet];
+ if ( calRoc ){
+ AliTRDCalROC* outlierROC = 0;
+ if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
+ if(type == 0) arr[n] = calRoc->GetRMS(outlierROC)*factor;
+ else arr[n] = calRoc->GetRMS(outlierROC);
+ n++;
+ }
+ }
+ return TMath::Mean(n,arr);
+}
+
+//_____________________________________________________________________________
+Double_t AliTRDCalPad::GetMedian(const AliTRDCalDet *calDet, Int_t type
+ , AliTRDCalPad* outlierPad)
+{
+ //
+ // return mean of the median of all ROCs
+ // If AliTRDCalDet, the correct the CalROC from the detector coefficient
+ // type == 0 for gain and vdrift
+ // type == 1 for t0
+ //
+ Double_t factor = 0.0;
+ if(type == 0) factor = 1.0;
+ Double_t arr[kNdet];
+ Int_t n=0;
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet) factor = calDet->GetValue(idet);
+ AliTRDCalROC *calRoc = fROC[idet];
+ if ( calRoc ){
+ AliTRDCalROC* outlierROC = 0;
+ if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
+ if(type == 0) arr[n] = calRoc->GetMedian(outlierROC)*factor;
+ else arr[n] = calRoc->GetMedian(outlierROC)+factor;
+ n++;
+ }
+ }
+ return TMath::Mean(n,arr);
+}
+
+//_____________________________________________________________________________
+Double_t AliTRDCalPad::GetLTM(Double_t *sigma, Double_t fraction
+ , const AliTRDCalDet *calDet, Int_t type
+ , AliTRDCalPad* outlierPad)
+{
+ //
+ // return mean of the LTM and sigma of all ROCs
+ // If calDet correct the CalROC from the detector coefficient
+ // type == 0 for gain and vdrift
+ // type == 1 for t0
+ //
+ Double_t factor = 0.0;
+ if(type == 0) factor = 1.0;
+ Double_t arrm[kNdet];
+ Double_t arrs[kNdet];
+ Double_t *sTemp=0x0;
+ Int_t n=0;
+
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet) factor = calDet->GetValue(idet);
+ AliTRDCalROC *calRoc = fROC[idet];
+ if ( calRoc ){
+ if ( sigma ) sTemp=arrs+n;
+ AliTRDCalROC* outlierROC = 0;
+ if (outlierPad) outlierROC = outlierPad->GetCalROC(idet);
+ if(type == 0) arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC)*factor;
+ else arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC)+factor;
+ n++;
+ }
+ }
+ if ( sigma ) *sigma = TMath::Mean(n,arrs);
+ return TMath::Mean(n,arrm);
+}
+
+//_____________________________________________________________________________
+TH1F * AliTRDCalPad::MakeHisto1D(const AliTRDCalDet *calDet, Int_t typedet
+ , Float_t min, Float_t max,Int_t type)
+{
+ //
+ // make 1D histo
+ // type -1 = user defined range
+ // 0 = nsigma cut nsigma=min
+ // If calDet correct the CalROC from the detector coefficient
+ // typedet == 0 for gain and vdrift
+ // typedet == 1 for t0
+ //
+ Float_t kEpsilonr = 0.005;
+
+ Double_t factor = 0.0;
+ if(typedet == 0) factor = 1.0;
+
+ if (type>=0){
+ if (type==0){
+ // nsigma range
+ Float_t mean = GetMean(calDet,typedet);
+ Float_t sigma = 0.0;
+ if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
+ else {
+ Double_t rms = 0.0;
+ sigma = GetMeanRMS(rms,calDet,typedet);
+ }
+ 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(calDet,typedet);
+ 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,calDet,typedet);
+ sigma*=min;
+ if(sigma < kEpsilonr) sigma = kEpsilonr;
+ min = mean-sigma;
+ max = mean+sigma;
+ }
+ }
+ char name[1000];
+ sprintf(name,"%s Pad 1D",GetTitle());
+ TH1F * his = new TH1F(name,name,100, min,max);
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet) factor = calDet->GetValue(idet);
+ if (fROC[idet]){
+ for (Int_t irow=0; irow<fROC[idet]->GetNrows(); irow++){
+ for (Int_t icol=0; icol<fROC[idet]->GetNcols(); icol++){
+ if(typedet == 0) his->Fill(fROC[idet]->GetValue(irow,icol)*factor);
+ else his->Fill(fROC[idet]->GetValue(irow,icol)+factor);
+ }
+ }
+ }
+ }
+ return his;
+}
+
+//_____________________________________________________________________________
+TH2F *AliTRDCalPad::MakeHisto2DSmPl(Int_t sm, Int_t pl, const AliTRDCalDet *calDet
+ , Int_t typedet, Float_t min, Float_t max,Int_t type)
+{
+ //
+ // Make 2D graph
+ // sm - supermodule number
+ // pl - plane number
+ // If calDet correct the CalROC from the detector coefficient
+ // typedet == 0 for gain and vdrift
+ // typedet == 1 for t0
+ //
+ gStyle->SetPalette(1);
+ Double_t factor = 0.0;
+ if(typedet == 0) factor = 1.0;
+
+ Float_t kEpsilon = 0.000000000001;
+ Float_t kEpsilonr = 0.005;
+
+ if (type>=0){
+ if (type==0){
+ // nsigma range
+ Float_t mean = GetMean(calDet,typedet);
+ Float_t sigma = 0.0;
+ if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
+ else {
+ Double_t rms = 0.0;
+ sigma = GetMeanRMS(rms,calDet,typedet);
+ }
+ 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(calDet,typedet);
+ 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,calDet,typedet);
+ sigma*=min;
+ if(sigma < kEpsilonr) sigma = kEpsilonr;
+ min = mean-sigma;
+ max = mean+sigma;
+ }
+ }
+
+ AliTRDpadPlane *padPlane0 = new AliTRDpadPlane(pl,0);
+ Double_t row0 = padPlane0->GetRow0();
+ Double_t col0 = padPlane0->GetCol0();
+
+ char name[1000];
+ sprintf(name,"%s Pad 2D sm %d pl %d",GetTitle(),sm,pl);
+ TH2F * his = new TH2F( name, name, 76,-TMath::Abs(row0),TMath::Abs(row0),144,-TMath::Abs(col0),TMath::Abs(col0));
+
+ // Where we begin
+ Int_t offsetsmpl = 30*sm+pl;
+
+
+ for (Int_t k = 0; k < kNcham; k++){
+ Int_t det = offsetsmpl+k*6;
+ if(calDet) factor = calDet->GetValue(det);
+ if (fROC[det]){
+ AliTRDCalROC * calRoc = fROC[det];
+ for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
+ for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
+ if (TMath::Abs(calRoc->GetValue(icol,irow))>kEpsilon){
+ Int_t binz = 0;
+ Int_t kb = kNcham-1-k;
+ if(kb > 2) binz = 16*(kb-1)+12+irow+1;
+ else binz = 16*kb+irow+1;
+ Int_t biny = icol+1;
+ Float_t value = calRoc->GetValue(icol,irow);
+ if(typedet == 0) his->SetBinContent(binz,biny,value*factor);
+ else his->SetBinContent(binz,biny,value+factor);
+ }
+ }
+ }
+ }
+ }
+ his->SetXTitle("z (cm)");
+ his->SetYTitle("y (cm)");
+ his->SetStats(0);
+ his->SetMaximum(max);
+ his->SetMinimum(min);
+ return his;
+}
+
+//_____________________________________________________________________________
+TH2F *AliTRDCalPad::MakeHisto2DCh(Int_t ch, const AliTRDCalDet *calDet, Int_t typedet
+ , Float_t min, Float_t max,Int_t type)
+{
+ //
+ // Make 2D graph mean value in z direction
+ // ch - chamber
+ // If calDet correct the CalROC from the detector coefficient
+ // typedet == 0 for gain and vdrift
+ // typedet == 1 for t0
+ //
+ gStyle->SetPalette(1);
+ Double_t factor = 0.0;
+ if(typedet == 0) factor = 1.0;
+
+ Float_t kEpsilonr = 0.005;
+
+ if (type>=0){
+ if (type==0){
+ // nsigma range
+ Float_t mean = GetMean(calDet,typedet);
+ Float_t sigma = 0.0;
+ if(GetRMS(calDet,typedet) > kEpsilonr) sigma = GetRMS(calDet,typedet);
+ else {
+ Double_t rms = 0.0;
+ sigma = GetMeanRMS(rms,calDet,typedet);
+ }
+ 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(calDet,typedet);
+ 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,calDet,typedet);
+ sigma*=min;
+ if(sigma < kEpsilonr) sigma = kEpsilonr;
+ min = mean-sigma;
+ max = mean+sigma;
+ }
+ }
+
+ AliTRDgeometry *TRDgeo = new AliTRDgeometry();
+ TRDgeo->Init();
+
+ Float_t kEpsilon = 0.000000000001;
+
+ Double_t poslocal[3] = {0.0,0.0,0.0};
+ Double_t posglobal[3] = {0.0,0.0,0.0};
+
+ char name[1000];
+ sprintf(name,"%s Pad 2D ch %d",GetTitle(),ch);
+ TH2F * his = new TH2F( name, name, 400,-400.0,400.0,400,-400.0,400.0);
+
+ // Where we begin
+ Int_t offsetch = 6*ch;
+
+
+ for (Int_t isec = 0; isec < kNsect; isec++){
+ for(Int_t ipl = 0; ipl < kNplan; ipl++){
+ Int_t det = offsetch+isec*30+ipl;
+ if(calDet) factor = calDet->GetValue(det);
+ if (fROC[det]){
+ AliTRDCalROC * calRoc = fROC[det];
+ AliTRDpadPlane *padPlane = new AliTRDpadPlane(ipl,ch);
+ for (Int_t icol=0; icol<calRoc->GetNcols(); icol++){
+ poslocal[0] = TRDgeo->GetTime0(ipl);
+ poslocal[2] = padPlane->GetRowPos(0);
+ poslocal[1] = padPlane->GetColPos(icol);
+ TRDgeo->RotateBack(det,poslocal,posglobal);
+ Int_t binx = 1+TMath::Nint((posglobal[0]+400.0)*0.5);
+ Int_t biny = 1+TMath::Nint((posglobal[1]+400.0)*0.5);
+ Float_t value = 0.0;
+ Int_t nb = 0;
+ for (Int_t irow=0; irow<calRoc->GetNrows(); irow++){
+ if (TMath::Abs(calRoc->GetValue(icol,irow))>kEpsilon){
+ value += calRoc->GetValue(icol,irow);
+ nb++;
+ }
+ }
+ value = value/nb;
+ if(typedet == 0) his->SetBinContent(binx,biny,value*factor);
+ else his->SetBinContent(binx,biny,value+factor);
+ }
+ }
+ }
+ }
+ his->SetXTitle("x (cm)");
+ his->SetYTitle("y (cm)");
+ his->SetStats(0);
+ his->SetMaximum(max);
+ his->SetMinimum(min);
+ return his;
+}
+
+//_____________________________________________________________________________
+Bool_t AliTRDCalPad::Add(Float_t c1)
+{
+ //
+ // add constant for all channels of all ROCs
+ //
+
+ Bool_t result = kTRUE;
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if (fROC[idet]){
+ if(!fROC[idet]->Add(c1)) result = kFALSE;
}
}
+ return result;
+}
+
+//_____________________________________________________________________________
+Bool_t AliTRDCalPad::Multiply(Float_t c1)
+{
+ //
+ // multiply constant for all channels of all ROCs
+ //
+ Bool_t result = kTRUE;
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if (fROC[idet]){
+ if(!fROC[idet]->Multiply(c1)) result = kFALSE;
+ }
+ }
+ return result;
+}
+//_____________________________________________________________________________
+Bool_t AliTRDCalPad::Add(const AliTRDCalPad * pad, Double_t c1
+ , const AliTRDCalDet * calDet1, const AliTRDCalDet *calDet2
+ , Int_t type)
+{
+ //
+ // add calpad channel by channel multiplied by c1 - all ROCs
+ // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
+ // then you have calDet1 and the calPad together
+ // calDet2 and pad together
+ // typedet == 0 for gain and vdrift
+ // typedet == 1 for t0
+ //
+ Float_t kEpsilon = 0.000000000001;
+
+ Double_t factor1 = 0.0;
+ Double_t factor2 = 0.0;
+ if(type == 0) {
+ factor1 = 1.0;
+ factor2 = 1.0;
+ }
+ Bool_t result = kTRUE;
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet1) factor1 = calDet1->GetValue(idet);
+ if(calDet2) factor2 = calDet2->GetValue(idet);
+ if (fROC[idet]){
+ if(type == 0){
+ if(TMath::Abs(factor1) > kEpsilon){
+ if(!fROC[idet]->Add(pad->GetCalROC(idet),c1*factor2/factor1)) result = kFALSE;
+ }
+ else result = kFALSE;
+ }
+ else{
+ AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
+ if(!croc->Add(factor2)) result = kFALSE;
+ if(!fROC[idet]->Add(croc,c1)) result = kFALSE;
+ }
+ }
+ }
+ return result;
}
+//_____________________________________________________________________________
+Bool_t AliTRDCalPad::Multiply(const AliTRDCalPad * pad, const AliTRDCalDet * calDet1
+ , const AliTRDCalDet *calDet2, Int_t type)
+{
+ //
+ // multiply calpad channel by channel - all ROCs
+ // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
+ // then you have calDet1 and the calPad together
+ // calDet2 and pad together
+ // typedet == 0 for gain and vdrift
+ // typedet == 1 for t0
+ //
+ Float_t kEpsilon = 0.000000000001;
+ Bool_t result = kTRUE;
+ Double_t factor1 = 0.0;
+ Double_t factor2 = 0.0;
+ if(type == 0){
+ factor1 = 1.0;
+ factor2 = 1.0;
+ }
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet1) factor1 = calDet1->GetValue(idet);
+ if(calDet2) factor2 = calDet2->GetValue(idet);
+ if (fROC[idet]){
+ if(type == 0){
+ if(TMath::Abs(factor1) > kEpsilon){
+ AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
+ if(!croc->Multiply(factor2)) result = kFALSE;
+ if(!fROC[idet]->Multiply(croc)) result = kFALSE;
+ }
+ else result = kFALSE;
+ }
+ else{
+ AliTRDCalROC *croc2 = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
+ if(!croc2->Add(factor2)) result = kFALSE;
+ if(!fROC[idet]->Add(factor1)) result = kFALSE;
+ if(!fROC[idet]->Multiply(croc2)) result = kFALSE;
+ if(!fROC[idet]->Add(-factor1)) result = kFALSE;
+ }
+ }
+ }
+ return result;
+}
+
+//_____________________________________________________________________________
+Bool_t AliTRDCalPad::Divide(const AliTRDCalPad * pad, const AliTRDCalDet * calDet1
+ , const AliTRDCalDet *calDet2, Int_t type)
+{
+ //
+ // divide calpad channel by channel - all ROCs
+ // If calDet1 and calDet2, the correct the CalROC from the detector coefficient
+ // then you have calDet1 and the calPad together
+ // calDet2 and pad together
+ // typedet == 0 for gain and vdrift
+ // typedet == 1 for t0
+ //
+ Float_t kEpsilon = 0.000000000001;
+ Bool_t result = kTRUE;
+ Double_t factor1 = 0.0;
+ Double_t factor2 = 0.0;
+ if(type == 0){
+ factor1 = 1.0;
+ factor2 = 1.0;
+ }
+ for (Int_t idet = 0; idet < kNdet; idet++) {
+ if(calDet1) factor1 = calDet1->GetValue(idet);
+ if(calDet2) factor2 = calDet2->GetValue(idet);
+ if (fROC[idet]){
+ if(type == 0){
+ if(TMath::Abs(factor1) > kEpsilon){
+ AliTRDCalROC *croc = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
+ if(!croc->Multiply(factor2)) result = kFALSE;
+ if(!fROC[idet]->Divide(croc)) result = kFALSE;
+ }
+ else result = kFALSE;
+ }
+ else{
+ AliTRDCalROC *croc2 = new AliTRDCalROC((const AliTRDCalROC) *pad->GetCalROC(idet));
+ if(!croc2->Add(factor2)) result = kFALSE;
+ if(!fROC[idet]->Add(factor1)) result = kFALSE;
+ if(!fROC[idet]->Divide(croc2)) result = kFALSE;
+ if(!fROC[idet]->Add(-factor1)) result = kFALSE;
+ }
+ }
+ }
+ return result;
+}
///////////////////////////////////////////////////////////////////////////////
// //
-// TRD calibration class for parameters which are saved per pad //
+// TRD calibration class for parameters which are saved per pad //
// //
///////////////////////////////////////////////////////////////////////////////
class AliTRDCalROC;
class AliTRDCalDet;
+class TH2F;
+class TH1F;
-class AliTRDCalPad : public TNamed {
+class AliTRDCalPad : public TNamed
+{
public:
AliTRDCalROC *GetCalROC(Int_t p, Int_t c, Int_t s) const
{ return fROC[GetDet(p,c,s)]; };
- void ScaleROCs(AliTRDCalDet* values);
+ Bool_t ScaleROCs(const AliTRDCalDet* values);
+
+ // Statistic
+ Double_t GetMeanRMS(Double_t &rms, const AliTRDCalDet *calDet = 0, Int_t type = 0);
+ Double_t GetMean(const AliTRDCalDet *calDet = 0, Int_t type = 0, AliTRDCalPad* outlierPad = 0);
+ Double_t GetRMS(const AliTRDCalDet *calDet = 0, Int_t type = 0, AliTRDCalPad* outlierPad = 0) ;
+ Double_t GetMedian(const AliTRDCalDet *calDet = 0, Int_t type = 0, AliTRDCalPad* outlierPad = 0) ;
+ Double_t GetLTM(Double_t *sigma=0, Double_t fraction=0.9, const AliTRDCalDet *calDet = 0, Int_t type = 0, AliTRDCalPad* outlierPad = 0);
+
+ // Plot functions
+ TH1F *MakeHisto1D(const AliTRDCalDet *calDet = 0, Int_t typedet=0, Float_t min=4, Float_t max=-4,Int_t type=0);
+ TH2F *MakeHisto2DSmPl(Int_t sm, Int_t pl, const AliTRDCalDet *calDet = 0, Int_t typedet=0, Float_t min=4, Float_t max=-4,Int_t type=0);
+ TH2F *MakeHisto2DCh(Int_t ch, const AliTRDCalDet *calDet = 0, Int_t typedet=0, Float_t min=4, Float_t max=-4,Int_t type=0);
+
+ // Algebra functions
+ Bool_t Add(Float_t c1);
+ Bool_t Multiply(Float_t c1);
+ Bool_t Add(const AliTRDCalPad * pad, Double_t c1 = 1, const AliTRDCalDet * calDet1 = 0, const AliTRDCalDet *calDet2 = 0, Int_t type = 0);
+ Bool_t Multiply(const AliTRDCalPad * pad, const AliTRDCalDet * calDet1 = 0, const AliTRDCalDet *calDet2 = 0, Int_t type = 0);
+ Bool_t Divide(const AliTRDCalPad * pad, const AliTRDCalDet * calDet1 = 0, const AliTRDCalDet *calDet2 = 0, Int_t type = 0);
protected:
// //
///////////////////////////////////////////////////////////////////////////////
+#include <iostream>
+#include <fstream>
+#include <string>
+#include <TStyle.h>
+
#include "AliTRDCalROC.h"
+#include "TMath.h"
+#include "AliMathBase.h"
+#include "TLinearFitter.h"
+#include "TArrayI.h"
+#include "TH2F.h"
+#include "TH1F.h"
+#include "TArrayF.h"
+#include "TGraph2D.h"
+#include "TGraphDelaunay.h"
+#include "TList.h"
+
+#include "AliTRDCommonParam.h"
+#include "AliTRDpadPlane.h"
+#include "AliLog.h"
ClassImp(AliTRDCalROC)
Int_t iBin = 0;
- if (((AliTRDCalROC &) c).fData) delete [] ((AliTRDCalROC &) c).fData;
- ((AliTRDCalROC &) c).fData = new UShort_t[fNchannels];
+ fData = new UShort_t[fNchannels];
for (iBin = 0; iBin < fNchannels; iBin++) {
- ((AliTRDCalROC &) c).fData[iBin] = fData[iBin];
+ fData[iBin] = ((AliTRDCalROC &) c).fData[iBin];
}
}
for (iBin = 0; iBin < fNchannels; iBin++) {
((AliTRDCalROC &) c).fData[iBin] = fData[iBin];
}
-
+
TObject::Copy(c);
}
-//_____________________________________________________________________________
-void AliTRDCalROC::Scale(Float_t value)
+//___________________________________________________________________________________
+Double_t AliTRDCalROC::GetMean(AliTRDCalROC* outlierROC)
+{
+ //
+ // 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)))) {
+ ddata[NPoints]= (Double_t) fData[NPoints]/10000;
+ NPoints++;
+ }
+ }
+ Double_t mean = TMath::Mean(NPoints,ddata);
+ delete [] ddata;
+ return mean;
+}
+
+//_______________________________________________________________________________________
+Double_t AliTRDCalROC::GetMedian(AliTRDCalROC* outlierROC)
+{
+ //
+ // 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)))) {
+ ddata[NPoints]= (Double_t) fData[NPoints]/10000;
+ NPoints++;
+ }
+ }
+ Double_t mean = TMath::Median(NPoints,ddata);
+ delete [] ddata;
+ return mean;
+}
+
+//____________________________________________________________________________________________
+Double_t AliTRDCalROC::GetRMS(AliTRDCalROC* outlierROC)
+{
+ //
+ // 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)))) {
+ ddata[NPoints]= (Double_t)fData[NPoints]/10000;
+ NPoints++;
+ }
+ }
+ Double_t mean = TMath::RMS(NPoints,ddata);
+ delete [] ddata;
+ return mean;
+}
+
+//______________________________________________________________________________________________
+Double_t AliTRDCalROC::GetLTM(Double_t *sigma, Double_t fraction, AliTRDCalROC* outlierROC)
{
//
- // Scales all values of this ROC with the provided parameter. Is used if ROC defines
- // local variations of a global (or per detector defined) parameter
+ // Calculate LTM mean and sigma
//
- for (Int_t iBin = 0; iBin < fNchannels; iBin++) {
- fData[iBin] = (UShort_t) (value * fData[iBin]);
+ 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))) {
+ ddata[NPoints]= (Double_t) fData[NPoints]/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;
+}
+
+//__________________________________________________________________________________
+TH2F * AliTRDCalROC::MakeHisto2D(Float_t min, Float_t max,Int_t type)
+{
+ //
+ // 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));
+ }
+ }
+ his->SetStats(0);
+ his->SetMaximum(max);
+ his->SetMinimum(min);
+ return his;
+}
+
+//_______________________________________________________________________________________
+TH1F * AliTRDCalROC::MakeHisto1D(Float_t min, Float_t max,Int_t type)
+{
+ //
+ // 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));
+ }
+ }
+ return his;
}
/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
* See cxx source for full Copyright notice */
-/* $Id: AliTRDCalROC.h,v */
+/* $Id$ */
//////////////////////////////////////////////////
// //
//////////////////////////////////////////////////
#include <TObject.h>
+#include <TMath.h>
+#include <TLinearFitter.h>
+
+class TArrayI;
+class TArrayF;
+class TH2F;
+class TH1F;
+class TGraph2D;
+class TH2D;
//_____________________________________________________________________________
-class AliTRDCalROC : public TObject {
+class AliTRDCalROC : public TObject
+{
public:
void SetValue(Int_t col, Int_t row, Float_t value)
{ SetValue(GetChannel(col,row), value); };
- void Scale(Float_t value);
-
+ // statistic
+ //
+ Double_t GetMean(AliTRDCalROC *outlierROC=0);
+ Double_t GetRMS(AliTRDCalROC *outlierROC=0);
+ Double_t GetMedian(AliTRDCalROC *outlierROC=0);
+ Double_t GetLTM(Double_t *sigma=0, Double_t fraction=0.9, AliTRDCalROC *outlierROC=0);
+
+ // algebra
+ Bool_t Add(Float_t c1);
+ Bool_t Multiply(Float_t c1);
+ Bool_t Add(const AliTRDCalROC * roc, Double_t c1 = 1);
+ Bool_t Multiply(const AliTRDCalROC * roc);
+ Bool_t Divide(const AliTRDCalROC * roc);
+
+ //Plots
+ TH2F * MakeHisto2D(Float_t min, Float_t max,Int_t type);
+ TH1F * MakeHisto1D(Float_t min, Float_t max,Int_t type);
+
protected:
Int_t fPla; // Plane number