ATO-78 - edge handling (mirroring out of range pads) set as an option
authormivanov <marian.ivanov@cern.ch>
Mon, 8 Sep 2014 08:16:14 +0000 (10:16 +0200)
committermivanov <marian.ivanov@cern.ch>
Tue, 28 Oct 2014 23:13:06 +0000 (00:13 +0100)
TPC/Base/AliTPCCalPad.cxx
TPC/Base/AliTPCCalPad.h
TPC/Base/AliTPCCalROC.cxx
TPC/Base/AliTPCCalROC.h

index 1b3ef41..df3c6c3 100644 (file)
@@ -159,27 +159,29 @@ void AliTPCCalPad::SetCalROC(AliTPCCalROC* roc, Int_t sector){
       fROC[sector]->SetValue(ichannel, roc->GetValue(ichannel));
 }
 
-Bool_t  AliTPCCalPad::MedianFilter(Int_t deltaRow, Int_t deltaPad){
+Bool_t  AliTPCCalPad::MedianFilter(Int_t deltaRow, Int_t deltaPad, AliTPCCalPad*outlierPad,  Bool_t doEdge){
   //
   // replace constent with median in the neigborhood
   //
   Bool_t isOK=kTRUE;
   for (Int_t isec = 0; isec < kNsec; isec++) {
+    AliTPCCalROC *outlierROC=(outlierPad==NULL)?NULL:outlierPad->GetCalROC(isec);
     if (fROC[isec]){
-      isOK&=fROC[isec]->MedianFilter(deltaRow,deltaPad);
+      isOK&=fROC[isec]->MedianFilter(deltaRow,deltaPad,outlierROC,doEdge);
     }
   }
   return isOK;
 }
 
-Bool_t  AliTPCCalPad::LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type){
+Bool_t  AliTPCCalPad::LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type, AliTPCCalPad*outlierPad,  Bool_t doEdge){
   //
   // replace constent with LTM statistic  in  neigborhood
   //
   Bool_t isOK=kTRUE;
   for (Int_t isec = 0; isec < kNsec; isec++) {
+    AliTPCCalROC *outlierROC=(outlierPad==NULL)?NULL:outlierPad->GetCalROC(isec);
     if (fROC[isec]){
-      isOK&=fROC[isec]->LTMFilter(deltaRow, deltaPad,fraction,type);
+      isOK&=fROC[isec]->LTMFilter(deltaRow, deltaPad,fraction,type,outlierROC,doEdge);
     }
   }
   return isOK;
@@ -913,6 +915,10 @@ AliTPCCalPad *AliTPCCalPad::MakePadFromTree(TTree * treePad, const char *query,
   //
   // make cal pad from the tree 
   //
+  if (!treePad){
+    ::Error("AliTPCCalPad::MakePadFromTree(TTree * treePad, const char *query, const char* name)","Input tree is missing");
+    return 0;
+  }
   if (treePad->GetEntries()!=kNsec) return 0;
   AliTPCCalPad * calPad= new AliTPCCalPad(name,name);
   if (name) calPad->SetName(name);
index 814718e..5898be2 100644 (file)
@@ -48,8 +48,8 @@ class AliTPCCalPad : public TNamed {
   void AddFriend(TTree * tree, const char *friendName, const char *fname=0);
   //
   // convolution
-  Bool_t MedianFilter(Int_t deltaRow, Int_t deltaPad);
-  Bool_t LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type);
+  Bool_t MedianFilter(Int_t deltaRow, Int_t deltaPad, AliTPCCalPad*outlierPad=0, Bool_t doEdge=kTRUE);
+  Bool_t LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type, AliTPCCalPad*outlierPad=0, Bool_t doEdge=kTRUE);
   //
   // algebra
   void Add(Float_t c1);   // add constant c1 to all channels of all ROCs
index 20ab069..ae3c991 100644 (file)
@@ -146,7 +146,8 @@ void AliTPCCalROC::Streamer(TBuffer &R__b)
 
 //
 
-Bool_t AliTPCCalROC::MedianFilter(Int_t deltaRow, Int_t deltaPad){
+Bool_t AliTPCCalROC::MedianFilter(Int_t deltaRow, Int_t deltaPad, AliTPCCalROC* outlierROC,  Bool_t doEdge){
+  //){
   //
   //   Modify content of the object - raplace value by median in neighorhood
   //
@@ -176,14 +177,19 @@ Bool_t AliTPCCalROC::MedianFilter(Int_t deltaRow, Int_t deltaPad){
          if (sign<0){
            kRow=iRow-dRow;
            jPad=iPad-dPad+offset;
+           if (!doEdge) continue;
          }       
          if (IsInRange(UInt_t(kRow),jPad)){
-           cacheBuffer[counter]=sign*(GetValue(kRow,jPad)-value);
-           counter++;
+           Bool_t isOutlier=(outlierROC==NULL)?kFALSE:outlierROC->GetValue(kRow,jPad)>0;
+           if (!isOutlier){
+             cacheBuffer[counter]=sign*(GetValue(kRow,jPad)-value);
+             counter++;
+           }
          }
        }
       }
-      newBuffer[fkIndexes[iRow]+iPad] = TMath::Median(counter, cacheBuffer)+value;
+      newBuffer[fkIndexes[iRow]+iPad]=0.;
+      if (counter>1) newBuffer[fkIndexes[iRow]+iPad] = TMath::Median(counter, cacheBuffer)+value;
     }
   }
   memcpy(fData, newBuffer,GetNchannels()*sizeof(Float_t));
@@ -194,7 +200,8 @@ Bool_t AliTPCCalROC::MedianFilter(Int_t deltaRow, Int_t deltaPad){
 
 
 
-Bool_t AliTPCCalROC::LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type){
+Bool_t AliTPCCalROC::LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type, AliTPCCalROC* outlierROC,  Bool_t doEdge){
+  //){
   //
   //
   // //
@@ -225,17 +232,21 @@ Bool_t AliTPCCalROC::LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction,
          if (jPad<0) sign=-1;
          if (jPad>=jPads) sign=-1;
          if (sign<0){
+           if (!doEdge) continue;
            kRow=iRow-dRow;
            jPad=iPad-dPad+offset;
          } 
          if (IsInRange(UInt_t(kRow),jPad)){
-           cacheBuffer[counter]=sign*(GetValue(kRow,jPad)-value);
-           counter++;
+           Bool_t isOutlier=(outlierROC==NULL)?kFALSE:outlierROC->GetValue(kRow,jPad)>0;
+           if (!isOutlier){
+             cacheBuffer[counter]=sign*(GetValue(kRow,jPad)-value);
+             counter++;
+           }
          }
        }
       }
-      Double_t mean,rms;
-      AliMathBase::EvaluateUni(counter,cacheBuffer,mean,rms,TMath::Nint(fraction*Double_t(counter)));
+      Double_t mean=0,rms=0;
+      if (TMath::Nint(fraction*Double_t(counter))>1 ) AliMathBase::EvaluateUni(counter,cacheBuffer,mean,rms,TMath::Nint(fraction*Double_t(counter)));
       mean+=value;
       newBuffer[fkIndexes[iRow]+iPad] = (type==0)? mean:rms;
     }
index 1baee19..a85a460 100644 (file)
@@ -39,8 +39,8 @@ class AliTPCCalROC : public TNamed {
   void         SetValue(UInt_t channel, Float_t vd) {fData[channel]= vd; };
   virtual void Draw(Option_t* option = "");
   //
-  Bool_t MedianFilter(Int_t deltaRow, Int_t deltaPad);
-  Bool_t LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type);
+  Bool_t MedianFilter(Int_t deltaRow, Int_t deltaPad, AliTPCCalROC*outlierROC=0, Bool_t doEdge=kTRUE);
+  Bool_t LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type,  AliTPCCalROC*outlierROC=0, Bool_t doEdge=kTRUE);
   //
   // algebra
   void Add(Float_t c1); // add c1 to each channel of the ROC