First version of the fast ALTRO decoding for TPC (Jens)
authorcvetan <cvetan@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 17 Sep 2007 13:53:51 +0000 (13:53 +0000)
committercvetan <cvetan@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 17 Sep 2007 13:53:51 +0000 (13:53 +0000)
TPC/AliTPCRawStreamFast.cxx [new file with mode: 0644]
TPC/AliTPCRawStreamFast.h [new file with mode: 0644]
TPC/TPCbaseLinkDef.h
TPC/libTPCbase.pkg
TPC/testRawReaderFastDDL.C [new file with mode: 0644]
TPC/testRawReaderFastTPC.C [new file with mode: 0644]

diff --git a/TPC/AliTPCRawStreamFast.cxx b/TPC/AliTPCRawStreamFast.cxx
new file mode 100644 (file)
index 0000000..ae7210c
--- /dev/null
@@ -0,0 +1,159 @@
+/**************************************************************************
+ * 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$ */
+
+///////////////////////////////////////////////////////////////////////////////
+///
+/// This class provides access to TPC digits in raw data.
+///
+/// It loops over all TPC digits in the raw data given by the AliRawReader.
+/// The Next method goes to the next digit. If there are no digits left
+/// it returns kFALSE.
+/// Several getters provide information about the current digit.
+///
+///////////////////////////////////////////////////////////////////////////////
+
+#include <TSystem.h>
+
+#include "AliTPCRawStreamFast.h"
+#include "AliRawReader.h"
+#include "AliLog.h"
+#include "AliTPCAltroMapping.h"
+
+ClassImp(AliTPCRawStreamFast)
+
+//_____________________________________________________________________________
+AliTPCRawStreamFast::AliTPCRawStreamFast(AliRawReader* rawReader) :
+  AliAltroRawStreamFast(rawReader),
+  fSector(-1),
+  fPrevSector(-1),
+  fRow(-1),
+  fPrevRow(-1),
+  fPad(-1),
+  fPrevPad(-1),
+  fIsMapOwner(kTRUE)
+{
+  // create an object to read TPC raw digits
+
+  SelectRawData("TPC");
+
+  TString path = gSystem->Getenv("ALICE_ROOT");
+  path += "/TPC/mapping/Patch";
+  TString path2;
+  for(Int_t i = 0; i < 6; i++) {
+    path2 = path;
+    path2 += i;
+    path2 += ".data";
+    fMapping[i] = new AliTPCAltroMapping(path2.Data());
+  }
+
+  //fNoAltroMapping = kFALSE;
+}
+/*
+//_____________________________________________________________________________
+AliTPCRawStreamFast::AliTPCRawStreamFast(const AliTPCRawStreamFast& stream) :
+  AliAltroRawStreamFast(stream),
+  fSector(stream.fSector),
+  fPrevSector(stream.fPrevSector),
+  fRow(stream.fRow),
+  fPrevRow(stream.fPrevRow),
+  fPad(stream.fPad),
+  fPrevPad(stream.fPrevPad),
+  fIsMapOwner(kFALSE)
+{
+  for(Int_t i = 0; i < 6; i++) fMapping[i] = stream.fMapping[i];
+}
+
+//_____________________________________________________________________________
+AliTPCRawStreamFast& AliTPCRawStreamFast::operator = (const AliTPCRawStreamFast& stream)
+{
+  if(&stream == this) return *this;
+
+  ((AliAltroRawStreamFast *)this)->operator=(stream);
+
+  fSector = stream.fSector;
+  fPrevSector = stream.fPrevSector;
+  fRow = stream.fRow;
+  fPrevRow = stream.fPrevRow;
+  fPad = stream.fPad;
+  fPrevPad = stream.fPrevPad;
+  fIsMapOwner = kFALSE;
+
+  for(Int_t i = 0; i < 6; i++) fMapping[i] = stream.fMapping[i];
+
+  return *this;
+}
+*/
+//_____________________________________________________________________________
+AliTPCRawStreamFast::~AliTPCRawStreamFast()
+{
+// destructor
+
+  if (fIsMapOwner)
+    for(Int_t i = 0; i < 6; i++) delete fMapping[i];
+}
+
+//_____________________________________________________________________________
+void AliTPCRawStreamFast::Reset()
+{
+  // reset tpc raw stream params
+  AliAltroRawStreamFast::Reset();
+  fSector = fPrevSector = fRow = fPrevRow = fPad = fPrevPad = -1;
+}
+
+//_____________________________________________________________________________
+Bool_t AliTPCRawStreamFast::NextChannel()
+{
+  // Read next TPC Channel
+  // Apply the TPC altro mapping to get
+  // the sector,pad-row and pad indeces
+  fPrevSector = fSector;
+  fPrevRow = fRow;
+  fPrevPad = fPad;
+  if (AliAltroRawStreamFast::NextChannel()) {
+      //    if (IsNewHWAddress())
+      if ( GetHWAddress() > -1 )
+      ApplyAltroMapping();
+    return kTRUE;
+  }
+  else
+    return kFALSE;
+}
+
+//_____________________________________________________________________________
+void AliTPCRawStreamFast::ApplyAltroMapping()
+{
+  // Take the DDL index, load
+  // the corresponding altro mapping
+  // object and fill the sector,row and pad indeces
+  Int_t ddlNumber = GetDDLNumber();
+  Int_t patchIndex;
+  if (ddlNumber < 72) {
+    fSector = ddlNumber / 2;
+    patchIndex = ddlNumber % 2;
+  }
+  else {
+    fSector = (ddlNumber - 72) / 4 + 36;
+    patchIndex = (ddlNumber - 72) % 4 + 2;
+  }
+
+  Short_t hwAddress = GetHWAddress();
+  fRow = fMapping[patchIndex]->GetPadRow(hwAddress);
+  fPad = fMapping[patchIndex]->GetPad(hwAddress);
+
+//  if ((fRow < 0) || (fPad < 0))
+//    AddMappingErrorLog(Form("hw=%d",hwAddress));
+}
diff --git a/TPC/AliTPCRawStreamFast.h b/TPC/AliTPCRawStreamFast.h
new file mode 100644 (file)
index 0000000..5b25d99
--- /dev/null
@@ -0,0 +1,55 @@
+#ifndef ALITPCRAWSTREAMFAST_H
+#define ALITPCRAWSTREAMFAST_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+///////////////////////////////////////////////////////////////////////////////
+///
+/// This class provides access to TPC digits in raw data.
+///
+///////////////////////////////////////////////////////////////////////////////
+
+#include "AliAltroRawStreamFast.h"
+
+class AliRawReader;
+class AliAltroMapping;
+
+class AliTPCRawStreamFast: public AliAltroRawStreamFast {
+  public :
+    AliTPCRawStreamFast(AliRawReader* rawReader);
+    virtual ~AliTPCRawStreamFast();
+
+    virtual void             Reset();
+    virtual Bool_t           NextChannel();
+
+    inline Int_t GetSector()     const { return fSector; }     // Provide index of current sector
+    inline Int_t GetPrevSector() const { return fPrevSector; } // Provide index of previous sector
+    inline Bool_t  IsNewSector() const {return fSector != fPrevSector;};
+    inline Int_t GetRow()        const { return fRow; }        // Provide index of current row
+    inline Int_t GetPrevRow()    const { return fPrevRow; }    // Provide index of previous row
+    inline Bool_t  IsNewRow()    const {return (fRow != fPrevRow) || IsNewSector();};
+    inline Int_t GetPad()        const { return fPad; }        // Provide index of current pad
+    inline Int_t GetPrevPad()    const { return fPrevPad; }    // Provide index of previous pad
+    inline Bool_t  IsNewPad()    const {return (fPad != fPrevPad) || IsNewRow();};
+
+
+  protected :
+    AliTPCRawStreamFast& operator = (const AliTPCRawStreamFast& stream);
+    AliTPCRawStreamFast(const AliTPCRawStreamFast& stream);
+
+    virtual void ApplyAltroMapping();
+
+    Int_t            fSector;       // index of current sector
+    Int_t            fPrevSector;   // index of previous sector
+    Int_t            fRow;          // index of current row
+    Int_t            fPrevRow;      // index of previous row
+    Int_t            fPad;          // index of current pad
+    Int_t            fPrevPad;      // index of previous pad
+
+    AliAltroMapping *fMapping[6];   // Pointers to ALTRO mapping
+    Bool_t           fIsMapOwner;   // does object own its mappings?
+
+    ClassDef(AliTPCRawStreamFast, 0)    // base class for reading TPC raw digits using the fast algorithm
+};
+
+#endif
index 194f9cc..5dcde9c 100644 (file)
@@ -44,6 +44,7 @@
 
 #pragma link C++ class AliTPCAltroMapping+;
 #pragma link C++ class AliTPCRawStream+;
+#pragma link C++ class AliTPCRawStreamFast+;
 #pragma link C++ class AliTPCCalibPedestal+;
 #pragma link C++ class AliTPCCalibPulser+;
 #pragma link C++ class AliTPCCalibCE+;
index 8b8b611..6558956 100644 (file)
@@ -7,7 +7,7 @@ SRCS:=  AliSegmentID.cxx  AliSegmentArray.cxx AliDigits.cxx AliH2F.cxx \
        AliTPCmapper.cxx \
        AliTPCROC.cxx AliTPCCalROC.cxx AliTPCCalPad.cxx AliTPCCalDet.cxx \
        AliTPCcalibDB.cxx \
-       AliTPCAltroMapping.cxx AliTPCRawStream.cxx \
+       AliTPCAltroMapping.cxx AliTPCRawStream.cxx AliTPCRawStreamFast.cxx \
        AliTPCLaserTracks.cxx AliTPCSensorTemp.cxx  AliTPCSensorTempArray.cxx \
        AliTPCCalibPedestal.cxx AliTPCCalibPulser.cxx   AliTPCCalibCE.cxx \
         AliTPCPreprocessor.cxx  \
diff --git a/TPC/testRawReaderFastDDL.C b/TPC/testRawReaderFastDDL.C
new file mode 100644 (file)
index 0000000..51bdd58
--- /dev/null
@@ -0,0 +1,144 @@
+#include <stdio.h>
+#include <TString.h>
+#include <TROOT.h>
+#include <TStopwatch.h>
+#include <TH2I.h>
+#include <TFile.h>
+
+#include "AliRawReader.h"
+#include "AliRawReaderRoot.h"
+#include "AliLog.h"
+
+#include "AliAltroRawStreamFast.h"
+#include "AliAltroRawStream.h"
+
+/*
+ compare old and new AltroRawStream algorithms for each pad and timebin
+
+ check if bins are filled twice (which should not be the case!!!)
+*/
+
+void testRawReaderFastDDL(const Char_t *file="/data.local/data/06000002142000.1A.root")
+{
+    // set minimal screen output
+    AliLog::SetGlobalDebugLevel(0) ;
+    AliLog::SetGlobalLogLevel(AliLog::kFatal);
+
+//    TString filename("/d/alice05/testtpc/raw/pulser/06000002142000.1A.root");
+    TString filename(file);
+
+    printf("File: %s\n", filename.Data());
+
+    AliRawReader *rawReader = new AliRawReaderRoot(filename);
+    if ( !rawReader ) return;
+    rawReader->RewindEvents();
+
+    AliAltroRawStreamFast *sf = new AliAltroRawStreamFast(rawReader);
+    AliAltroRawStream      *s = new AliAltroRawStream(rawReader);
+
+    s->SetOldRCUFormat(kTRUE);
+    s->SetNoAltroMapping(kFALSE);
+    s->SelectRawData("TPC");
+
+    Int_t ievent = 0;
+    Int_t count=0;
+
+    TH2I *h2ddlT1[216];
+    TH2I *h2ddlT2[216];
+
+    for ( Int_t i=0; i<216; i++ ){
+       h2ddlT1[i] = 0x0;
+       h2ddlT2[i] = 0x0;
+    }
+
+    TStopwatch timer1;
+    TStopwatch timer2;
+    TStopwatch timer3;
+
+    while (rawReader->NextEvent()){
+       printf("\nevent: %d\n",ievent);
+       Bool_t input=kFALSE;
+
+
+       //old algorithm
+       timer1.Start();timer2.Start(kFALSE);
+       while ( s->Next() ){
+            if ( !h2ddlT1[s->GetDDLNumber()] ) h2ddlT1[s->GetDDLNumber()] = new TH2I(Form("hddl1_%d",s->GetDDLNumber()),"h2c1",3584,0,3584,1024,0,1024);
+           TH2I *hist = h2ddlT1[s->GetDDLNumber()];
+
+           //fast filling, TH1::Fill takes awfully long
+           Int_t bin=(s->GetTime()+1)*(3584+2)+s->GetHWAddress()+1;
+           // check if this bin was allready filled
+           if ( hist->GetArray()[bin] > 0 )
+               printf(" not 0:   |  %.3d : %.3d (%.3d)\n",
+                       s->GetHWAddress(), s->GetSignal(), hist->GetArray()[bin]);
+           else
+               hist->GetArray()[bin]=s->GetSignal();
+
+
+           input=kTRUE;
+           count++;
+       }
+       timer1.Stop();timer2.Stop();
+       printf("old  --  Time: %.4f (%.4f)\n", timer1.RealTime(), timer1.CpuTime());
+       // END OLD
+
+       rawReader->Reset();
+
+       //new algorithm
+        timer1.Start();timer3.Start(kFALSE);
+       while ( sf->NextDDL() ){
+           if ( !h2ddlT2[sf->GetDDLNumber()] ) h2ddlT2[sf->GetDDLNumber()] = new TH2I(Form("hddl2_%d",s->GetDDLNumber()),"h2c1",3584,0,3584,1024,0,1024);
+           TH2I *hist = h2ddlT2[sf->GetDDLNumber()];
+           while ( sf->NextChannel() ){
+               UInt_t signal=0;
+               while ( sf->NextBunch() ){
+                   for (UInt_t timebin=sf->GetStartTimeBin(); timebin<sf->GetEndTimeBin(); timebin++){
+                       signal=sf->GetSignals()[timebin-sf->GetStartTimeBin()];
+
+                       //fast filling, TH1::Fill takes awfully long
+                       Int_t bin=(timebin+1+1)*(3584+2)+sf->GetHWAddress()+1; // timebins of old and new algorithm differ by 1!!!
+                        // check if this bin was allready filled
+                       if ( hist->GetArray()[bin] > 0 )
+                           printf(" not 0:   |  %.3d : %.3d (%.3d)\n",
+                                   sf->GetHWAddress(), signal, hist->GetArray()[bin]);
+                       else
+                           hist->GetArray()[bin]=signal;
+                   }
+               }
+           }
+       }
+       timer1.Stop();timer3.Stop();
+       printf("new  --  Time: %.4f (%.4f)\n", timer1.RealTime(), timer1.CpuTime());
+       // END NEW
+
+        //check if all data are the same for both algorithms
+       for ( Int_t ddl=0; ddl<216; ddl++ ){
+           TH2I *hist = h2ddlT1[ddl];
+           if ( !hist ) continue;
+           TH2I *hist2 = h2ddlT2[ddl];
+           for ( Int_t hadd=0; hadd<3584; hadd++ )
+               for ( Int_t time=0; time<1024; time++ ){
+                   Int_t bin=(time+1)*(3584+2)+hadd+1;
+                   Int_t val1 = hist->GetArray()[bin];
+                   Int_t val2 = hist2->GetArray()[bin];
+                   if ( val1 != val2 )
+                       printf("%.2d. %.3d %.4d %.4d: %d - %d = %d\n", ievent, ddl, hadd, time, val1, val2, val1-val2);
+
+                    //reset for the next event
+                   hist->GetArray()[bin]=0;
+                   hist2->GetArray()[bin]=0;
+               }
+       }
+
+       if (input) ievent++;
+    }
+
+    printf("total old  --  Time: %.4f (%.4f)\n", timer2.RealTime(), timer2.CpuTime());
+    printf("total new  --  Time: %.4f (%.4f)\n", timer3.RealTime(), timer3.CpuTime());
+
+    delete rawReader;
+    delete [] h2ddlT1;
+    delete [] h2ddlT2;
+}
+
diff --git a/TPC/testRawReaderFastTPC.C b/TPC/testRawReaderFastTPC.C
new file mode 100644 (file)
index 0000000..b668bb4
--- /dev/null
@@ -0,0 +1,207 @@
+#include <stdio.h>
+#include <TString.h>
+#include <TROOT.h>
+#include <TStopwatch.h>
+#include <TCanvas.h>
+#include <TFile.h>
+#include <TH2F.h>
+
+#include "AliLog.h"
+
+#include "AliRawReader.h"
+#include "AliRawReaderRoot.h"
+
+#include "AliTPCRawStream.h"
+#include "AliTPCRawStreamFast.h"
+
+#include "AliAltroRawStream.h"
+#include "AliAltroRawStreamFast.h"
+
+#include "AliTPCCalPad.h"
+#include "AliTPCCalROC.h"
+
+/*
+    .L testRawReaderFastTPC.C+
+    testRawReaderFastTPC()
+*/
+
+
+AliTPCCalPad *testRawReaderFastTPC(const Char_t *file="/data.local/data/06000002142000.1A.root")
+{
+
+    AliLog::SetGlobalDebugLevel(0) ;
+    AliLog::SetGlobalLogLevel(AliLog::kFatal);
+//    TString filename("/d/alice05/testtpc/raw/pulser/06000002142000.1A.root");  //nfs
+//    TString filename("root://lxfs35.gsi.de:1094//alice/testtpc/raw2006/06000001537001.001.root");
+    TString filename(file);  //local
+    // on castor: /castor/cern.ch/alice/data/2006/09/18/15/06000002142000.1A.root
+
+
+    printf("File: %s\n", filename.Data());
+
+    AliRawReader *rawReader = new AliRawReaderRoot(filename);
+    if ( !rawReader ) return 0x0;
+    rawReader->RewindEvents();
+
+    AliTPCRawStreamFast *sf = new AliTPCRawStreamFast(rawReader);
+    AliTPCRawStream      *s = new AliTPCRawStream(rawReader);
+
+    s->SetOldRCUFormat(kTRUE);
+    s->SetNoAltroMapping(kFALSE);
+
+    Int_t ievent = 0;
+    Int_t ievent2 = 0;
+    Int_t count=0;
+
+    AliTPCCalPad *padOld = new AliTPCCalPad("old","old");
+    AliTPCCalPad *padNew = new AliTPCCalPad("new","new");
+    AliTPCCalPad *padSum = new AliTPCCalPad("sum","sum");
+    AliTPCCalPad *padhadd = new AliTPCCalPad("sum","sum");
+
+    for (Int_t sec=0; sec<72; sec++){
+       AliTPCCalROC *roc = padhadd->GetCalROC(sec);
+       for (UInt_t ch=0; ch<roc->GetNchannels(); ch++)
+           roc->SetValue(ch,-1);
+    }
+
+    TStopwatch timer1;
+    TStopwatch timer2;
+    TStopwatch timer3;
+
+    TCanvas *c1 = (TCanvas*)gROOT->FindObject("c1");
+    if ( !c1 ) c1 = new TCanvas("c1","c1");
+    c1->Clear();
+    c1->Divide(2,2);
+
+
+    while (rawReader->NextEvent()){
+       //old algorithm
+       printf("\nevent: %d (%d)\n",ievent, ievent2);
+
+       Bool_t input=kFALSE;
+
+//     if ( ievent != 19 ) input=kTRUE;
+//     else{
+           timer1.Start();timer2.Start(kFALSE);
+
+           Int_t sum1 = 0;
+            Int_t sum2 = 0;
+       while ( s->Next() ){
+            AliTPCCalROC *roc = padhadd->GetCalROC(s->GetSector());
+
+            //check if hwaddress gets overwritten
+           if ( roc->GetValue(s->GetRow(), s->GetPad()) == -1 )
+               roc->SetValue(s->GetRow(), s->GetPad(), s->GetHWAddress());
+           else
+               if ( roc->GetValue(s->GetRow(), s->GetPad()) != s->GetHWAddress()){
+                   printf("#%.2d: %.2d.%.3d.%.2d,%.4d - old [%.4d] now[%.4d]",
+                          ievent, s->GetSector(), s->GetRow(),s->GetPad(), s->GetTime(),
+                          (Int_t)roc->GetValue(s->GetRow(), s->GetPad()), s->GetHWAddress());
+               }
+
+
+           Float_t val=padOld->GetCalROC(s->GetSector())->GetValue(s->GetRow(), s->GetPad());
+           padOld->GetCalROC(s->GetSector())->SetValue(s->GetRow(), s->GetPad(), s->GetSignal()+val);
+
+/*         if ( ievent == 19 && s->GetSector() == 25 && s->GetRow() == 00 && s->GetPad()==9 ){
+               printf("old         | (%.2d.%.3d.%.2d.%.4d | %.3d ): %.3d\n",
+                      s->GetSector(), s->GetRow(), s->GetPad(), s->GetTime(), s->GetHWAddress(), s->GetSignal());
+                sum1+=s->GetSignal();
+           }
+*/
+
+           input=kTRUE;
+           count++;
+       }
+       timer1.Stop();timer2.Stop();
+       printf("old  --  Time: %.4f (%.4f)\n", timer1.RealTime(), timer1.CpuTime());
+
+       rawReader->Reset();
+
+       //new algorithm
+        timer1.Start();timer3.Start(kFALSE);
+       while ( sf->NextDDL() ){
+           while ( sf->NextChannel() ){
+               UInt_t signal=0;
+               while ( sf->NextBunch() ){
+                   for (UInt_t timebin=sf->GetStartTimeBin(); timebin<sf->GetEndTimeBin(); timebin++){
+
+                       AliTPCCalROC *roc = padhadd->GetCalROC(sf->GetSector());
+//                     if ( roc->GetValue(sf->GetRow(), sf->GetPad()) == -1 )
+//                         roc->SetValue(sf->GetRow(), sf->GetPad(), sf->GetHWAddress());
+//                     else
+                           if ( roc->GetValue(sf->GetRow(), sf->GetPad()) != sf->GetHWAddress()){
+                               printf("#%.2d: %.2d.%.3d.%.2d,%.4d - old [%.4d] now[%.4d]",
+                                      ievent, sf->GetSector(), sf->GetRow(),sf->GetPad(), timebin+1,
+                                      (Int_t)roc->GetValue(sf->GetRow(), sf->GetPad()), sf->GetHWAddress());
+                           }
+                            /*
+                           Int_t sig = sf->GetSignals()[timebin-sf->GetStartTimeBin()];
+                           if ( ievent == 19 && sf->GetSector() == 25 && sf->GetRow() == 00 && sf->GetPad()==9  ){
+                               printf("new         | (%.2d.%.3d.%.2d.%.4d | %.3d ): %.3d\n",
+                                      sf->GetSector(), sf->GetRow(), sf->GetPad(), timebin+1, sf->GetHWAddress(), sig );
+                               sum2+=sig;
+                           }
+                            */
+                       signal+=sf->GetSignals()[timebin-sf->GetStartTimeBin()];
+                   }
+                   padNew->GetCalROC(sf->GetSector())->SetValue(sf->GetRow(),sf->GetPad(),signal);
+               }
+           }
+       }
+       timer1.Stop();timer3.Stop();
+       printf("new  --  Time: %.4f (%.4f)\n", timer1.RealTime(), timer1.CpuTime());
+
+       printf("sum1: %d, sum2: %d, diff: %d\n",sum1,sum2,sum1-sum2);
+
+       AliTPCCalPad com(*padOld);
+       com.Add(padNew, -1);
+
+       c1->cd(1);
+       padOld->MakeHisto2D(1)->Draw("colz");
+       c1->cd(2);
+       padNew->MakeHisto2D(1)->Draw("colz");
+       c1->cd(3);
+       com.MakeHisto2D(1)->Draw("colz");
+
+        //loop over all sectors, rows, pads
+       for ( UInt_t iSec=0; iSec<72; iSec++ )
+           for ( UInt_t iRow=0; iRow<com.GetCalROC(iSec)->GetNrows(); iRow++ )
+               for ( UInt_t iPad=0; iPad<com.GetCalROC(iSec)->GetNPads(iRow); iPad++){
+                   Float_t val = com.GetCalROC(iSec)->GetValue(iRow, iPad);
+                   Float_t valo = padNew->GetCalROC(iSec)->GetValue(iRow, iPad);
+                    Float_t valn = padOld->GetCalROC(iSec)->GetValue(iRow, iPad);
+                    //check if values for old and new algorithm differ
+//                 if ( val != 0 ){
+                   if ( (Int_t)valo != (Int_t)valn ){
+                        Float_t hadd = padhadd->GetCalROC(iSec)->GetValue(iRow,iPad);
+                       printf("Event: %.2d | (%.2d.%.3d.%.2d | %f ): %.3f\n", ievent, iSec, iRow, iPad, hadd, val);
+                       padSum->GetCalROC(iSec)->SetValue(iRow, iPad, padSum->GetCalROC(iSec)->GetValue(iRow, iPad)+1);
+                   }
+                   com.GetCalROC(iSec)->SetValue(iRow, iPad, 0);
+                   padOld->GetCalROC(iSec)->SetValue(iRow, iPad, 0);
+                   padNew->GetCalROC(iSec)->SetValue(iRow, iPad, 0);
+               }
+
+
+       c1->cd(4);
+       padSum->MakeHisto2D(1)->Draw("colz");
+
+       c1->Modified();
+       c1->Update();
+  //      }// end event sel
+       if (input) ievent++;
+       ievent2++;
+    }
+    TFile f("output.root","recreate");
+    padSum->Write();
+    f.Save();
+    f.Close();
+
+    printf("total old  --  Time: %.4f (%.4f)\n", timer2.RealTime(), timer2.CpuTime());
+    printf("total new  --  Time: %.4f (%.4f)\n", timer3.RealTime(), timer3.CpuTime());
+
+    delete rawReader;
+    return padSum;
+}
+