Merged Cvetans RowHoughTransformer, Anders latest developments in comp
authorloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 30 Mar 2004 18:34:26 +0000 (18:34 +0000)
committerloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 30 Mar 2004 18:34:26 +0000 (18:34 +0000)
and src, as well as little changes concering the build process.

104 files changed:
HLT/Makefile
HLT/Makefile.conf
HLT/Makefile.rules
HLT/bin/sethlt.csh
HLT/bin/sethlt.sh
HLT/bin/sethlt_cern.csh
HLT/bin/sethlt_cern.sh
HLT/comp/AliL3ClusterFitter.cxx
HLT/comp/AliL3ClusterFitter.h
HLT/comp/AliL3Compress.cxx
HLT/comp/AliL3DataCompressor.cxx
HLT/comp/AliL3DataCompressorHelper.cxx
HLT/comp/AliL3ModelTrack.cxx
HLT/comp/AliL3Modeller.cxx
HLT/comp/AliL3OfflineDataCompressor.cxx
HLT/doc/README
HLT/doc/changelog.bin
HLT/doc/changelog.comp
HLT/doc/changelog.doc
HLT/doc/changelog.exa
HLT/doc/changelog.hough
HLT/doc/changelog.kalman
HLT/doc/changelog.misc
HLT/doc/changelog.programs
HLT/doc/changelog.src
HLT/doc/changelog.top
HLT/doc/changelog.trigger
HLT/doc/taginfo
HLT/exa/evalrowhough.C [new file with mode: 0644]
HLT/exa/evaltracker.C
HLT/exa/rootlogon.C
HLT/exa/runrowhough.C [new file with mode: 0644]
HLT/exa/runtracker.C
HLT/hough/AliL3Histogram.cxx
HLT/hough/AliL3Histogram.h
HLT/hough/AliL3Histogram1D.cxx
HLT/hough/AliL3HistogramAdaptive.cxx
HLT/hough/AliL3Hough.cxx
HLT/hough/AliL3Hough.h
HLT/hough/AliL3HoughBaseTransformer.cxx
HLT/hough/AliL3HoughBaseTransformer.h
HLT/hough/AliL3HoughClusterTransformer.cxx
HLT/hough/AliL3HoughDisplay.cxx
HLT/hough/AliL3HoughEval.cxx
HLT/hough/AliL3HoughEval.h
HLT/hough/AliL3HoughIntMerger.cxx
HLT/hough/AliL3HoughLinkDef.h
HLT/hough/AliL3HoughMaxFinder.cxx
HLT/hough/AliL3HoughMaxFinder.h
HLT/hough/AliL3HoughMerger.cxx
HLT/hough/AliL3HoughTest.cxx
HLT/hough/AliL3HoughTrack.cxx
HLT/hough/AliL3HoughTransformer.cxx
HLT/hough/AliL3HoughTransformerGap.cxx [deleted file]
HLT/hough/AliL3HoughTransformerGlobal.cxx
HLT/hough/AliL3HoughTransformerLUT.cxx
HLT/hough/AliL3HoughTransformerNew.cxx
HLT/hough/AliL3HoughTransformerRow.cxx [new file with mode: 0644]
HLT/hough/AliL3HoughTransformerRow.h [moved from HLT/hough/AliL3HoughTransformerGap.h with 54% similarity]
HLT/hough/AliL3HoughTransformerVhdl.cxx
HLT/hough/Makefile
HLT/kalman/AliL3Kalman.cxx
HLT/kalman/AliL3Kalman.h
HLT/kalman/AliL3KalmanTrack.cxx
HLT/kalman/AliL3KalmanTrack.h
HLT/misc/AliL3DDLDataFileHandler.cxx
HLT/misc/AliL3DDLRawReaderFile.cxx
HLT/misc/AliL3DataHandler.cxx
HLT/misc/AliL3Stopwatch.cxx
HLT/misc/AliL3TPCMapping.cxx
HLT/misc/AliL3TransBit.cxx
HLT/misc/AliL3VHDLClusterFinder.cxx
HLT/programs/Makefile
HLT/programs/ali2raw.cxx
HLT/programs/convcosmicsfile.cxx
HLT/programs/gettransform.cxx
HLT/programs/read.cxx
HLT/programs/runhough.cxx
HLT/programs/runit.cxx
HLT/programs/runtracker.cxx
HLT/programs/runvhdlcf.cxx
HLT/programs/runvhdlhough.cxx
HLT/src/AliL3ClustFinderNew.cxx
HLT/src/AliL3ConfMapFit.cxx
HLT/src/AliL3ConfMapTrack.cxx
HLT/src/AliL3Display.cxx
HLT/src/AliL3Evaluate.cxx
HLT/src/AliL3FileHandler.cxx
HLT/src/AliL3Fitter.cxx
HLT/src/AliL3Logger.h
HLT/src/AliL3Logging.h
HLT/src/AliL3RawDataFileHandler.cxx
HLT/src/AliL3RootTypes.h
HLT/src/AliL3StandardIncludes.h
HLT/src/AliL3Track.cxx
HLT/src/AliL3TrackArray.cxx
HLT/src/AliL3TrackSegmentData.h
HLT/src/AliL3Transform.cxx
HLT/src/AliLevel3.cxx
HLT/src/AliLevel3.h
HLT/trigger/AliD0Trigger.cxx
HLT/trigger/AliD0Trigger.h
HLT/trigger/RunD0Trigger.C
HLT/trigger/RunD0offline.C [new file with mode: 0644]

index 5e4bc0f..26f4094 100644 (file)
@@ -41,6 +41,7 @@ print:
        @echo "ALIHLT_NOLOGGING  = $(ALIHLT_NOLOGGING)"
        @echo "ALIHLT_DOMC       = $(ALIHLT_DOMC)"
        @echo "ALIHLT_ALIDETECT  = $(ALIHLT_ALIDETECT)"
+       @echo "ALIHLT_ROWHOUGH   = $(ALIHLT_ROWHOUGH)"
 
 help:
        cat doc/README
index 16b2913..7ca26ed 100644 (file)
@@ -54,25 +54,21 @@ ALIHLT_USENEWIO = false
 ifeq ($(ALICE_LEVEL),ali-head)
 ALIHLT_USENEWIO = true
 endif
-ifeq ($(ALICE_LEVEL),ali-v4-01-00)
-ALIHLT_USENEWIO = true
 endif
+
+ifeq ($(ALIHLT_ROWHOUGH),true)
+USEROWHOUGH = 1
 endif
 
 ifeq ($(ALIHLT_USENEWIO),true)
 USENEWIO = 1
 endif
 
-
 #----------------------------------------------------
 #Some compiler flags or defines: You can use your
 #own setting by defining them outside (make -e)
 #----------------------------------------------------
 
-GCCVERSION   = $(shell $(CXX) --version | head -n 1 | cut -d" " -f 3 | cut -d. -f 1 | cut -d" " -f1)
-CXXGCC3FLAGS = -DGCCVERSION=$(GCCVERSION)
-PROFILEFLAGS = -g -pg
-
 ifeq ($(ARCH),Darwin)
 FINKDIR         = /sw
 CXX             = g++
@@ -86,14 +82,15 @@ SOFLAGS         = -bundle -flat_namespace -undefined suppress
 DYFLAGS         = -dynamiclib -flat_namespace -undefined suppress \
                   -compatibility_version 1 -current_version 1.0.0 
 LDFLAGS         = -O $(EXTRALDFLAGS) -L/sw/lib -ldl
-LDSTATIC     = ar
-STATICFLAGS  = rucs 
 else
 CXX          = g++
 CXXFLAGS     = -O2 -fPIC -Wall $(CXXGCC3FLAGS) $(EXTRACXXFLAGS) 
 LD           = $(CXX)
 LDFLAGS      = -O2 $(EXTRALDFLAGS)
 SOFLAGS      = -shared
-LDSTATIC     = ar
-STATICFLAGS  = rucs
 endif
+
+#static flags for profiling
+PROFILEFLAGS = -g -pg
+LDSTATIC     = ar
+STATICFLAGS  = rucs 
index acbd4df..d201eef 100644 (file)
@@ -45,6 +45,10 @@ ifeq ($(DOMC),1)
 DEFSTR += -Ddo_mc
 endif
 
+ifeq ($(USEROWHOUGH),1)
+DEFSTR += -DROWHOUGH
+endif
+
 ifneq ($(NOLOGGING),1)
 DEFSTR += -Duse_logging
 ifdef ALIHLT_MLUCDIR
@@ -96,7 +100,7 @@ $(ALIHLT_STATIC): $(STATICOBJS)
 
 $(DICT): $(HDRS)
        @echo "Generating dictionary..."
-       rootcint -f $(DICT) -c $(CINTCXXFLAGS) $(INCLUDES) -DGCCVERSION=$(GCCVERSION) \
+       rootcint -f $(DICT) -c $(CINTCXXFLAGS) $(INCLUDES) \
                     $(DEFSTR) -include AliL3StandardIncludes.h $(HDRS)
 
 $(OBJDIR)/%.o: %.cxx 
@@ -124,6 +128,7 @@ print:
        @echo "ALIHLT_NOLOGGING  = $(ALIHLT_NOLOGGING)"
        @echo "ALIHLT_DOMC       = $(ALIHLT_DOMC)"
        @echo "ALIHLT_ALIDETECT  = $(ALIHLT_ALIDETECT)"
+       @echo "ALIHLT_ROWHOUGH   = $(ALIHLT_ROWHOUGH)"
        @echo "ROOTSTR           = $(ROOTSTR)"
        @echo "ALIROOTSTR        = $(ALIROOTSTR)"
 
index 7373982..18f9d60 100755 (executable)
@@ -16,7 +16,7 @@ setenv ALIHLT_LIBDIR ${ALIHLT_TOPDIR}/lib_${ALIHLT_USEPACKAGE}
 setenv ALIHLT_NOLOGGING false
 setenv ALIHLT_DOMC true
 setenv ALIHLT_ALIDETECT true
-
+setenv ALIHLT_ROWHOUGH false
 setenv ALIHLT_MLUCDIR ${ALIHLT_BASEDIR}/kip/MLUC
 
 setenv ALIHLT_DATADIR /data1/head
index bea7588..c8c50a9 100755 (executable)
@@ -20,14 +20,13 @@ export ALIHLT_LIBDIR=$ALIHLT_TOPDIR/lib_$ALIHLT_USEPACKAGE
 export ALIHLT_NOLOGGING=false
 export ALIHLT_DOMC=true
 export ALIHLT_ALIDETECT=true
-
+export ALIHLT_ROWHOUGH=false
 export ALIHLT_MLUCDIR=/usr/local/kip/MLUC
 
 #export ALIHLT_DATADIR=/mnt/local/alidata/head
 #export ALIHLT_TRANSFORMFILE=$ALIHLT_DATADIR/l3transform.config
 #export ALIHLT_GEOPATH=$ALIDATADIR
 
-
 if test -z "$LD_LIBRARY_PATH"; then
   export LD_LIBRARY_PATH=$ALIHLT_MLUCDIR/lib:$ALIHLT_LIBDIR
 elif test -z "`echo $LD_LIBRARY_PATH | grep $ALIHLT_MLUCDIR/lib`"; 
index 813f17a..829ec49 100755 (executable)
@@ -12,11 +12,12 @@ setenv ALIHLT_LIBDIR ${ALIHLT_TOPDIR}/lib_${ALIHLT_USEPACKAGE}
 
 setenv ALIHLT_NOLOGGING false
 setenv ALIHLT_DOMC true
-setenv ALIHLT_ALIDETECT true
+setenv ALIHLT_ALIDETECT false
+setenv ALIHLT_ROWHOUGH true
 
 setenv ALIHLT_MLUCDIR ${ALIHLT_BASEDIR}/kip/MLUC
 
-setenv ALIHLT_DATADIR /data1/head
+#setenv ALIHLT_DATADIR /data1/head
 #setenv ALIHLT_TRANSFORMFILE ${ALIHLT_DATADIR}/l3transform.config
 #setenv ALIHLT_GEOPATH ${ALIHLT_DATADIR}
 
index b2333e4..f5b3a77 100755 (executable)
@@ -13,7 +13,7 @@ export ALIHLT_LIBDIR=$ALIHLT_TOPDIR/lib_$ALIHLT_USEPACKAGE
 export ALIHLT_NOLOGGING=true
 export ALIHLT_DOMC=true
 export ALIHLT_ALIDETECT=false
-
+export ALIHLT_ROWHOUGH=true
 #export ALIHLT_MLUCDIR=/usr/local/kip/MLUC
 
 #export ALIHLT_DATADIR=/mnt/local/alidata/head
index b4002a9..e6aba41 100644 (file)
 #include "AliL3SpacePointData.h"
 #include "AliL3Compress.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
+/** /class AliL3ClusterFitter
+//<pre>
 //_____________________________________________________________
 //
 //  AliL3ClusterFitter
 //
-
+</pre>
+*/
 
 ClassImp(AliL3ClusterFitter)
 
@@ -217,10 +220,9 @@ void AliL3ClusterFitter::LoadLocalSegments()
     }
 }
 
-void AliL3ClusterFitter::LoadSeeds(Int_t *rowrange,Bool_t offline,Int_t eventnr)
+void AliL3ClusterFitter::LoadSeeds(Int_t *rowrange,Bool_t offline,Int_t eventnr,Float_t zvertex)
 {
   //Function assumes _global_ tracks written to a single file.
-  
   cout<<"Loading the seeds"<<endl;
   Char_t fname[1024];
   fEvent = eventnr;
@@ -293,6 +295,7 @@ void AliL3ClusterFitter::LoadSeeds(Int_t *rowrange,Bool_t offline,Int_t eventnr)
              //exit(5);
            }
          Float_t xyz_cross[3] = {track->GetPointX(),track->GetPointY(),track->GetPointZ()};
+         xyz_cross[2] += zvertex;
          
          Int_t sector,row;
          AliL3Transform::Slice2Sector(slice,j,sector,row);
@@ -340,6 +343,7 @@ void AliL3ClusterFitter::LoadSeeds(Int_t *rowrange,Bool_t offline,Int_t eventnr)
              xyz_cross[0] = track->GetPointX();
              xyz_cross[1] = track->GetPointY();
              xyz_cross[2] = track->GetPointZ();
+             xyz_cross[2] += zvertex;
              Int_t sector,row;
              AliL3Transform::Slice2Sector(slice,j,sector,row);
              AliL3Transform::Global2Raw(xyz_cross,sector,row);
@@ -1248,4 +1252,3 @@ void AliL3ClusterFitter::WriteClusters(Bool_t global)
   fClusters=0;
   fNClusters=0;
 }
-
index 1248c6b..33f42c8 100644 (file)
@@ -50,7 +50,7 @@ class AliL3ClusterFitter : public AliL3Modeller {
   
   void Init(Int_t slice,Int_t patch,Int_t *rowrange,AliL3TrackArray *tracks);
   void Init(Int_t slice,Int_t patch);
-  void LoadSeeds(Int_t *rowrange,Bool_t offline=kTRUE,Int_t eventnr=0);
+  void LoadSeeds(Int_t *rowrange,Bool_t offline=kTRUE,Int_t eventnr=0,Float_t zvertex=0.0);
   void LoadLocalSegments();
   void FindClusters();
   void AddClusters();
@@ -66,7 +66,6 @@ class AliL3ClusterFitter : public AliL3Modeller {
   Float_t GetZWidthFactor() {return fCurrentPadRow < AliL3Transform::GetLastRow(1) ? fZInnerWidthFactor : fZOuterWidthFactor;}
   AliL3TrackArray *GetSeeds() {return fSeeds;}
   
-  
   ClassDef(AliL3ClusterFitter,1) 
 
 };
index aff2b40..6df6cf4 100644 (file)
@@ -31,7 +31,7 @@
 
 #include "AliL3Compress.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 9325ba2..778a409 100644 (file)
@@ -40,7 +40,7 @@
 #include "AliL3DataCompressorHelper.h"
 #include "AliL3DataCompressor.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 3717041..4c2e34d 100644 (file)
@@ -10,7 +10,7 @@
 
 #include "AliL3DataCompressorHelper.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 212a6e6..dced528 100644 (file)
@@ -12,7 +12,7 @@
 
 #include "AliL3ModelTrack.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index f4553ab..29a1f5d 100644 (file)
@@ -18,7 +18,7 @@
 #include "AliL3FileHandler.h"
 #endif
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 2d2c4ea..2d3f52a 100644 (file)
@@ -33,7 +33,7 @@
 
 #include "AliL3OfflineDataCompressor.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 722515f..dd8cde7 100644 (file)
@@ -40,8 +40,9 @@ $ALIHLT_MLUCDIR --> directory of MLUC include and lib sub-directories,
 
 Compiliation flags:
 $ALIHLT_USENEWIO (false)  --> switch to NEWIO
-$ALIHLT_NOLOGGING (false) --> use HLT logger classes (requires MLUC lib)
+$ALIHLT_ROWHOUGH (false)  --> use row hough transformer
 $ALIHLT_DOMC (true)       --> store Monte Carlo ids (for ALIROOT package)
+$ALIHLT_NOLOGGING (false) --> use HLT logger classes (requires MLUC lib)
 $ALIHLT_ALIDETECT (true)  --> detect ALIROOT version using cvs
 
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
index 30c5175..cd07fb0 100644 (file)
@@ -1,3 +1,8 @@
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/bin/sethlt_cern.csh, /alice/cvs/hltcvs/level3code/bin/sethlt_cern.sh, /alice/cvs/hltcvs/level3code/bin/sethlt.csh, /alice/cvs/hltcvs/level3code/bin/sethlt.sh:
+       Added ALIHLT_ROWHOUGH variable.
+
 2003-07-29  Constantin Loizides  <loizides@hansa00>
 
        * /alice/cvs/hltcvs/level3code/bin/sethlt_cern.csh, /alice/cvs/hltcvs/level3code/bin/sethlt_cern.sh, /alice/cvs/hltcvs/level3code/bin/sethlt.csh, /alice/cvs/hltcvs/level3code/bin/sethlt.sh, /alice/cvs/hltcvs/level3code/bin/usehlt.csh, /alice/cvs/hltcvs/level3code/bin/usehlt.sh:
index bc8e3d8..f895837 100644 (file)
@@ -1,3 +1,12 @@
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/comp/AliL3ClusterFitter.cxx, /alice/cvs/hltcvs/level3code/comp/AliL3Compress.cxx, /alice/cvs/hltcvs/level3code/comp/AliL3DataCompressor.cxx, /alice/cvs/hltcvs/level3code/comp/AliL3DataCompressorHelper.cxx, /alice/cvs/hltcvs/level3code/comp/AliL3Modeller.cxx, /alice/cvs/hltcvs/level3code/comp/AliL3ModelTrack.cxx, /alice/cvs/hltcvs/level3code/comp/AliL3OfflineDataCompressor.cxx:
+       Removed GCCVERSION string from code. Replaced by compiler internal
+       macro __GNUC__ which is 3 for gcc version >3.
+
+       * /alice/cvs/hltcvs/level3code/comp/AliL3ClusterFitter.cxx, /alice/cvs/hltcvs/level3code/comp/AliL3ClusterFitter.h:
+       Added zvertex.
+
 2004-02-25  Anders Strand Vestbo  <vestbo@hansa00>
 
        * /alice/cvs/hltcvs/level3code/comp/AliL3Compress.cxx, /alice/cvs/hltcvs/level3code/comp/AliL3Compress.h:
index 29a07ed..7d755ab 100644 (file)
@@ -1,3 +1,10 @@
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/doc/README: Added ALIHLT_ROWHOUGH info.
+
+       * /alice/cvs/hltcvs/level3code/doc/taginfo:
+       Added new tag before big Cvetan changes.
+
 2004-01-22  Constantin Loizides  <loizides@hansa00>
 
        * /alice/cvs/hltcvs/level3code/doc/taginfo: Added new tag.
index cb75858..fc87092 100644 (file)
@@ -1,3 +1,31 @@
+2004-03-30  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/exa/evaltracker.C: Make it compilable.
+
+       * /alice/cvs/hltcvs/level3code/exa/runrowhough.C:
+       Little changes to have default params.
+
+       * /alice/cvs/hltcvs/level3code/exa/runtracker.C:
+       Save benchmark with bname.
+
+       * /alice/cvs/hltcvs/level3code/exa/evalrowhough.C:
+       Tried to make it compilable, but failled.
+
+2004-03-28  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/exa/runtracker.C:
+       Added name for tracker eval.
+
+       * /alice/cvs/hltcvs/level3code/exa/runrowhough.C:
+       Changes to compile macro.
+
+       * /alice/cvs/hltcvs/level3code/exa/rootlogon.C: Added runrowhough.
+
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/exa/evalrowhough.C, /alice/cvs/hltcvs/level3code/exa/runrowhough.C:
+       Added examples for new row hough transform.
+
 2004-02-02  Constantin Loizides  <loizides@hansa00>
 
        * /alice/cvs/hltcvs/level3code/exa/deconvclusters.C, /alice/cvs/hltcvs/level3code/exa/rootlogon.C, /alice/cvs/hltcvs/level3code/exa/runhough.C, /alice/cvs/hltcvs/level3code/exa/runtracker.C, /alice/cvs/hltcvs/level3code/exa/SetFitParameters.C, /alice/cvs/hltcvs/level3code/exa/SetHoughParameters.C:
index 485991b..b8b2d69 100644 (file)
@@ -1,3 +1,48 @@
+2004-03-28  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3HoughBaseTransformer.h:
+       Added old type for label.
+
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3Histogram1D.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HistogramAdaptive.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3Histogram.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughClusterTransformer.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3Hough.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughDisplay.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughEval.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughIntMerger.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughMaxFinder.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughMerger.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTest.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTrack.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformer.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerGlobal.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerLUT.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerNew.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerRow.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerVhdl.cxx:
+       Removed GCCVERSION string from code. Replaced by compiler internal
+       macro __GNUC__ which is 3 for gcc version >3.
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3HoughEval.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughEval.h:
+       Added zvertex.
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3Hough.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3Hough.h:
+       Added new row transformer (using special peakfinder method)
+       Added in addition support for reading data directly from
+       DATE, so that this version of hough transform can run directly
+       on the GDCs or during some data challenges for testing.
+       (changed methods init and constructor)
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3Histogram.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3Histogram.h:
+       Added GetPreciseBinCenter functions.
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3HoughBaseTransformer.cxx:
+       Cosmetics.
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3HoughMaxFinder.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughMaxFinder.h:
+       Added FindAdaptedRowPeaks function for row transformer.
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3HoughTrack.cxx:
+       Added ROWHOUGH def for mc label.
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3HoughBaseTransformer.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughBaseTransformer.h:
+       Added zvertex information.
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3HoughLinkDef.h, /alice/cvs/hltcvs/level3code/hough/Makefile:
+       Added new row transformer and removed old gap transformer.
+
+       * /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerGap.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerGap.h, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerRow.cxx, /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerRow.h:
+       Added Cvetans new version of fast Hough Transform called "Counting
+       Gaps and Rows." The idea is to count the number of consecutive rows
+       over gaps per entry of the Hough space. (The old version called
+       AliL3HoughTransformGap is removed.)
+
 2004-02-12  Constantin Loizides  <loizides@hansa00>
 
        * /alice/cvs/hltcvs/level3code/hough/AliL3HoughTransformerGap.cxx:
index 5e42f76..d1f1fc8 100644 (file)
@@ -1,3 +1,11 @@
+2004-03-18  Thomas Vik  <tvik@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/kalman/AliL3KalmanTrack.cxx:
+       Change in calling AliL3Transfrom::Global2Local command.
+
+       * /alice/cvs/hltcvs/level3code/kalman/AliL3KalmanTrack.h, /alice/cvs/hltcvs/level3code/kalman/AliL3KalmanTrack.cxx, /alice/cvs/hltcvs/level3code/kalman/AliL3Kalman.h, /alice/cvs/hltcvs/level3code/kalman/AliL3Kalman.cxx:
+       Changed global coordinates of padrows to local. Removed obsolete MakeSeed function.
+
 2004-02-02  Thomas Vik  <tvik@hansa00>
 
        * /alice/cvs/hltcvs/level3code/kalman/AliL3Kalman.cxx:
index 57426c1..e9b0a8f 100644 (file)
@@ -1,3 +1,12 @@
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/misc/AliL3DataHandler.cxx, /alice/cvs/hltcvs/level3code/misc/AliL3DDLDataFileHandler.cxx, /alice/cvs/hltcvs/level3code/misc/AliL3DDLRawReaderFile.cxx, /alice/cvs/hltcvs/level3code/misc/AliL3Stopwatch.cxx, /alice/cvs/hltcvs/level3code/misc/AliL3TPCMapping.cxx, /alice/cvs/hltcvs/level3code/misc/AliL3TransBit.cxx, /alice/cvs/hltcvs/level3code/misc/AliL3VHDLClusterFinder.cxx:
+       Removed GCCVERSION string from code. Replaced by compiler internal
+       macro __GNUC__ which is 3 for gcc version >3.
+
+       * /alice/cvs/hltcvs/level3code/misc/AliL3DDLRawReaderFile.cxx:
+       Added namespace for gcc-3.x
+
 2004-02-03  Constantin Loizides  <loizides@hansa00>
 
        * /alice/cvs/hltcvs/level3code/misc/AliL3DDLRawReaderFile.cxx:
index 2907e94..0ba3e60 100644 (file)
@@ -1,3 +1,9 @@
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/programs/ali2raw.cxx, /alice/cvs/hltcvs/level3code/programs/convcosmicsfile.cxx, /alice/cvs/hltcvs/level3code/programs/gettransform.cxx, /alice/cvs/hltcvs/level3code/programs/Makefile, /alice/cvs/hltcvs/level3code/programs/read.cxx, /alice/cvs/hltcvs/level3code/programs/runhough.cxx, /alice/cvs/hltcvs/level3code/programs/runit.cxx, /alice/cvs/hltcvs/level3code/programs/runtracker.cxx, /alice/cvs/hltcvs/level3code/programs/runvhdlcf.cxx, /alice/cvs/hltcvs/level3code/programs/runvhdlhough.cxx:
+       Removed GCCVERSION string from code. Replaced by compiler internal
+       macro __GNUC__ which is 3 for gcc version >3.
+
 2004-02-25  Constantin Loizides  <loizides@hansa00>
 
        * /alice/cvs/hltcvs/level3code/programs/gettransform.cxx, /alice/cvs/hltcvs/level3code/programs/Makefile, /alice/cvs/hltcvs/level3code/programs/runhough.cxx:
index 1707abe..31e7a47 100644 (file)
@@ -1,3 +1,43 @@
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/src/AliL3Transform.cxx:
+       Update version information.
+
+       * /alice/cvs/hltcvs/level3code/src/AliL3ClustFinderNew.cxx, /alice/cvs/hltcvs/level3code/src/AliL3Display.cxx, /alice/cvs/hltcvs/level3code/src/AliL3Evaluate.cxx, /alice/cvs/hltcvs/level3code/src/AliL3FileHandler.cxx, /alice/cvs/hltcvs/level3code/src/AliL3Logger.h, /alice/cvs/hltcvs/level3code/src/AliL3Logging.h, /alice/cvs/hltcvs/level3code/src/AliL3RawDataFileHandler.cxx, /alice/cvs/hltcvs/level3code/src/AliL3RootTypes.h, /alice/cvs/hltcvs/level3code/src/AliL3StandardIncludes.h, /alice/cvs/hltcvs/level3code/src/AliL3Track.cxx, /alice/cvs/hltcvs/level3code/src/AliL3Transform.cxx:
+       Removed GCCVERSION string from code. Replaced by compiler internal
+       macro __GNUC__ which is 3 for gcc version >3.
+
+       * /alice/cvs/hltcvs/level3code/src/AliL3TrackArray.cxx, /alice/cvs/hltcvs/level3code/src/AliL3TrackSegmentData.h:
+       Added ROWHOUGH flag to store weight and track id
+       in the track segment data structure. This is used
+       for the new row transformer and probably only
+       a temporary solution, as the proper way would be
+       to have a new AliL3HoughTrackRow class.
+
+       By default the flag $ALIHLT_ROWHOUGH points to false,
+       so there is no change to earlier versions.
+
+2004-03-15  Anders Strand Vestbo  <vestbo@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/src/AliL3ConfMapFit.cxx, /alice/cvs/hltcvs/level3code/src/AliL3ConfMapTrack.cxx, /alice/cvs/hltcvs/level3code/src/AliL3Fitter.cxx, /alice/cvs/hltcvs/level3code/src/AliLevel3.cxx, /alice/cvs/hltcvs/level3code/src/AliLevel3.h:
+       Bugfix related to track fit parameters <-> global track merging.
+
+       It turned out that there was a slight problem with the global track
+       merging in the case a vertex constraint was imposed in the circle fit.
+       The problem is quite involved, and is due to inconsistency when storing
+       the socalled first point on the track in AliL3Track::*fFirstPoint. This
+       point is set to the first point lying on the fit, but in the case of vertex
+       constraint this point is correspondingly set to the point of closest
+       approach to the vertex. HOWEVER, the global merger always assumes that
+       this point is the first associated cluster on the track, and thus we have
+       a conflict.
+
+       The (temporary) solution to the problem is to set the *fFirstPoint to the
+       innermost cluster on the track ALWAYS, i.e. in both vertexconstraint and
+       no vertex constraint. This means also that the linear fit in (s,z) space
+       does not currently include the vertex in the fit, in order to also set
+       Z0 to the innermost cluster of the track.
+
 2004-02-12  Constantin Loizides  <loizides@hansa00>
 
        * /alice/cvs/hltcvs/level3code/src/AliL3Transform.cxx, /alice/cvs/hltcvs/level3code/src/AliL3Transform.h:
index 9586833..f5ab958 100644 (file)
@@ -1,3 +1,19 @@
+2004-03-20  Constantin Loizides  <loizides@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/Makefile.conf, /alice/cvs/hltcvs/level3code/Makefile.rules:
+       Removed GCCVERSION string from code. Replaced by compiler internal
+       macro __GNUC__ which is 3 for gcc version >3.
+
+       * /alice/cvs/hltcvs/level3code/Makefile.conf, /alice/cvs/hltcvs/level3code/Makefile.rules:
+       Added ROWHOUGH flag to store weight and track id
+       in the track segment data structure. This is used
+       for the new row transformer and probably only
+       a temporary solution, as the proper way would be
+       to have a new AliL3HoughTrackRow class.
+
+       By default the flag $ALIHLT_ROWHOUGH points to false,
+       so there is no change to earlier versions.
+
 2004-02-25  Constantin Loizides  <loizides@hansa00>
 
        * /alice/cvs/hltcvs/level3code/Makefile.conf, /alice/cvs/hltcvs/level3code/Makefile.rules:
index 010f7af..436cad3 100644 (file)
@@ -1,3 +1,13 @@
+2004-03-12  Gaute Ovrebekk  <ovrebekk@hansa00>
+
+       * /alice/cvs/hltcvs/level3code/trigger/RunD0Trigger.C: small chages
+
+       * /alice/cvs/hltcvs/level3code/trigger/AliD0Trigger.cxx, /alice/cvs/hltcvs/level3code/trigger/AliD0Trigger.h:
+       Added more methodes
+
+       * /alice/cvs/hltcvs/level3code/trigger/RunD0offline.C:
+       Running with offline methids
+
 2004-02-24  Gaute Ovrebekk  <ovrebekk@hansa00>
 
        * /alice/cvs/hltcvs/level3code/trigger/AliD0Trigger.cxx:
index eeaa74d..abf97ba 100644 (file)
@@ -5,3 +5,7 @@ v1-1: AliTransform changes because of  Cosmics
 v1-2: First real release to CERN (with AliRoot 3-09-Release)
 v1-3: Deconvoluter and Hough changes tested (with AliRoot 3-09-Release)
 v1-4: Merge with cern-hlt tree to support NEWIO
+v1-5: Fitter and Track updates and status before intro of MACOSX handling
+v1-6: Compression and tracking updates as of March, 20. 
+
+
diff --git a/HLT/exa/evalrowhough.C b/HLT/exa/evalrowhough.C
new file mode 100644 (file)
index 0000000..2c09cce
--- /dev/null
@@ -0,0 +1,400 @@
+// $Id$
+
+struct GoodTrack 
+{
+  Int_t label;
+  Double_t eta;
+  Int_t code;
+  Double_t px,py,pz;
+  Double_t x,y,z;
+  Int_t nhits;
+  Int_t sector;
+};
+typedef struct GoodTrack GoodTrack;
+
+//Histograms
+TNtuple *fNtuppel;
+TH1F *fPtRes;
+TH1F *fNGoodTracksPt;
+TH1F *fNFoundTracksPt;
+TH1F *fNFakeTracksPt;
+TH1F *fTrackEffPt;
+TH1F *fFakeTrackEffPt;
+TH1F *fNGoodTracksEta;
+TH1F *fNFoundTracksEta;
+TH1F *fNFakeTracksEta;
+TH1F *fTrackEffEta;
+TH1F *fFakeTrackEffEta;
+TH1F *fNGoodTracksPhi;
+TH1F *fNFoundTracksPhi;
+TH1F *fNFakeTracksPhi;
+TH1F *fTrackEffPhi;
+TH1F *fFakeTrackEffPhi;
+TProfile *fFakesPt;
+TProfile *fFakesPhi;
+TProfile *fSlicesPt;
+TProfile *fSlicesPhi;
+TProfile2D *fFakesPtvsPhi;
+TNtuple *fNtupleRes;
+TH1F *fPtRes2;
+TH1F *fEtaRes;
+TH1F *fPsiRes;
+
+TH1F *fNGoodTracksRad;
+TH1F *fNFoundTracksRad;
+TH1F *fNGoodTracksZ;
+TH1F *fNFoundTracksZ;
+TH1F *fRadEff;
+TH1F *fZEff;
+
+void evalrowhough(Char_t *path="./fitter",Char_t *offlinepath="./",Int_t ievent=0, Float_t ptmin=0.4)  
+{
+  CreateHistos(25,0.1,4.1);
+
+  Char_t fname[1024];
+  sprintf(fname,"%s/tracks_%d.raw",path,ievent);
+  AliL3FileHandler *tfile = new AliL3FileHandler();
+  if(!tfile->SetBinaryInput(fname)){
+    LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open")
+      <<"Inputfile "<<fname<<" does not exist"<<ENDLOG; 
+    return;
+  }
+
+  AliL3TrackArray *fTracks = new AliL3TrackArray("AliL3HoughTrack");
+  tfile->Binary2TrackArray(fTracks);
+  //fTracks->QSort();
+  tfile->CloseBinaryInput();
+  delete tfile;
+
+  for(Int_t i=0; i<fTracks->GetNTracks(); i++)
+    {
+      AliL3HoughTrack *track = fTracks->GetCheckedTrack(i);
+      if(!track) continue; 
+      
+      track->SetEta(fTracks->GetCheckedTrack(i)->GetPseudoRapidity());
+      UInt_t *ids = track->GetHitNumbers();
+      Int_t slice = (ids[0]>>25)&0x7f;
+      track->SetSlice(slice);
+      cout<<"Track with label "<<ids<<" "<<track->GetNHits()<<" "<<track->GetMCid()<<" "<<ids[0]<<" "<<track->GetWeight()<<" "<<track->GetSlice()<<endl;
+    }
+
+  Char_t filename[1024];
+  sprintf(filename,"%s/good_tracks_tpc_%d",offlinepath,ievent);
+  ifstream in(filename);
+  if(!in)
+    {
+      cerr<<"AliL3Evaluate::GetGoodParticles : Problems opening file :"<<filename<<endl;
+      return;
+    }
+
+  Int_t MaxTracks=20000;
+  Int_t fGoodGen=0;
+  fGoodTracks = new GoodTrack[MaxTracks];
+  
+  while (in>>fGoodTracks[fGoodGen].label>>fGoodTracks[fGoodGen].code>>
+        fGoodTracks[fGoodGen].px>>fGoodTracks[fGoodGen].py>>fGoodTracks[fGoodGen].pz>>
+        fGoodTracks[fGoodGen].x>>fGoodTracks[fGoodGen].y>>fGoodTracks[fGoodGen].z)
+    {
+      fGoodTracks[fGoodGen].nhits=-1;
+      fGoodTracks[fGoodGen].sector=-1; 
+      fGoodGen++;
+      if (fGoodGen==MaxTracks) 
+       {
+         cerr<<"AliL3Evaluate::GetGoodParticles : Too many good tracks !\n";
+         break;
+       }
+    }
+
+  Int_t fGoodGen2 = 0;
+  Int_t fGoodFound = 0;
+  Int_t fGoodFound2 = 0;
+  Float_t fMinGoodPt = ptmin;
+  Float_t fMaxGoodPt = 99999.;
+
+  if(!fGoodTracks)
+    {
+      cerr<<"AliL3Evaluate::FillEffHistos : No good tracks"<<endl;
+      return;
+    }
+  cout<<"Comparing "<<fGoodGen<<" good tracks ..."<<endl;
+  for(Int_t i=0; i<fGoodGen; i++)
+    {
+      //cout<<"Checking particle "<<i<<endl;
+      Float_t ptg = TMath::Sqrt(fGoodTracks[i].px*fGoodTracks[i].px + fGoodTracks[i].py*fGoodTracks[i].py);
+      Float_t pg = TMath::Sqrt(fGoodTracks[i].px*fGoodTracks[i].px + fGoodTracks[i].py*fGoodTracks[i].py + fGoodTracks[i].pz*fGoodTracks[i].pz);
+      if(ptg < fMinGoodPt || ptg > fMaxGoodPt) continue;
+      Float_t pzg=fGoodTracks[i].pz;
+      Float_t dipangle=TMath::ATan2(pzg,ptg)*180./TMath::Pi();
+      Float_t etag=0.5 * TMath::Log((pg+pzg)/(pg-pzg));
+      Float_t phig=TMath::ATan2(fGoodTracks[i].py,fGoodTracks[i].px);
+      Float_t rg = TMath::Sqrt(fGoodTracks[i].x*fGoodTracks[i].x + fGoodTracks[i].y*fGoodTracks[i].y);
+      Float_t zg = TMath::Abs(fGoodTracks[i].z);
+      if (phig<0) phig+=2*TMath::Pi();
+
+      fNGoodTracksPt->Fill(ptg);
+      fNGoodTracksEta->Fill(etag);
+      fNGoodTracksPhi->Fill(phig);
+      if((zg != 0) || (rg != 0)) { 
+       fNGoodTracksRad->Fill(rg);
+       fNGoodTracksZ->Fill(zg);
+      }
+      fGoodGen2++;
+      Int_t found = 0;
+      Int_t nslices = 1;
+      Int_t org_slice;
+      for(Int_t k=0; k<fTracks->GetNTracks(); k++)
+       {
+         AliL3HoughTrack *track = (AliL3HoughTrack *)fTracks->GetCheckedTrack(k);
+         if(!track) continue;
+         //      Int_t nHits = track->GetNumberOfPoints();
+         //      if(nHits < fMinPointsOnTrack) break;
+         Int_t nHits = 1;
+         Int_t tracklabel;
+         tracklabel = track->GetMCid();
+         
+         if(TMath::Abs(tracklabel) != fGoodTracks[i].label) continue;
+         if(found == 0) {
+           found=1;
+           Float_t pt=track->GetPt();
+           org_slice = track->GetSlice();
+           Float_t eta = track->GetEta();
+           Float_t phi = track->GetPsi();
+           //cout<<eta<<" "<<phi<<" "<<etag<<" "<<phig<<" "<<track->GetEtaIndex()<<endl;
+           if(tracklabel == fGoodTracks[i].label) 
+             {
+               fGoodFound++;
+               
+               fNFoundTracksPt->Fill(ptg); 
+               fNFoundTracksEta->Fill(etag);
+               fNFoundTracksPhi->Fill(phig);
+               if((zg != 0) || (rg != 0)) { 
+                 fNFoundTracksRad->Fill(rg);
+                 fNFoundTracksZ->Fill(zg);
+               }
+               fNtuppel->Fill(ptg,pt,nHits);
+               fPtRes->Fill((pt-ptg)/ptg*100.);
+               fEtaRes->Fill(eta-etag);
+               fPsiRes->Fill(phi-phig);
+             }
+           else 
+             {
+               fNFakeTracksPt->Fill(ptg); 
+               fNFakeTracksEta->Fill(etag);
+               fNFakeTracksPhi->Fill(phig);
+             }
+           //fPtRes->Fill((pt-ptg)/ptg*100.);
+           //fNtuppel->Fill(ptg,pt,nHits);
+         }
+         else {
+           found++;
+           if(track->GetSlice() != org_slice)
+             nslices = 2;
+         }
+       }
+      if(!found)
+       //cout<<"Track "<<fGoodTracks[i].label<<" was not found"<<endl;
+       continue;
+      else {
+       fFakesPt->Fill(ptg,(Double_t)found-1.0);
+       fFakesPhi->Fill(phig,(Double_t)found-1.0);
+       fSlicesPt->Fill(ptg,(Double_t)nslices);
+       fSlicesPhi->Fill(phig,(Double_t)nslices);
+       //Float_t local_phi = phi-(phi/)
+       //fFakesPtvsPhi->Fill(1.0/ptg,
+      }
+      fGoodFound2 += found;
+    }
+
+  cout<<fGoodFound<<"("<<fGoodFound2<<") tracks found out of "<<fGoodGen2<<" generated tracks"<<endl;
+  cout<<" The overall efficiency is "<<(Float_t)fGoodFound/(Float_t)fGoodGen2<<endl;
+  
+  CalcEffHistos();
+  plotptres(ptmin);
+  Write2File();
+}
+
+void CreateHistos(Int_t nbin,Float_t xlow,Float_t xup)
+{
+  //Create the histograms 
+  
+  fNtuppel = new TNtuple("fNtuppel","Pt resolution","pt_gen:pt_found:nHits");
+  fNtuppel->SetDirectory(0);
+  fPtRes = new TH1F("fPtRes","Relative Pt resolution",30,-10.,10.); 
+  fNGoodTracksPt = new TH1F("fNGoodTracksPt","Good tracks vs pt",nbin,xlow,xup);    
+  fNFoundTracksPt = new TH1F("fNFoundTracksPt","Found tracks vs pt",nbin,xlow,xup);
+  fNFakeTracksPt = new TH1F("fNFakeTracksPt","Fake tracks vs pt",nbin,xlow,xup);
+  fTrackEffPt = new TH1F("fTrackEffPt","Tracking efficiency vs pt",nbin,xlow,xup);
+  fFakeTrackEffPt = new TH1F("fFakeTrackEffPt","Efficiency for fake tracks vs pt",nbin,xlow,xup);
+  
+  fNGoodTracksEta = new TH1F("fNGoodTracksEta","Good tracks vs eta",nbin,-1,1);
+  fNFoundTracksEta = new TH1F("fNFoundTracksEta","Found tracks vs eta",nbin,-1,1);
+  fNFakeTracksEta = new TH1F("fNFakeTracksEta","Fake tracks vs eta",nbin,-1,1);
+  fTrackEffEta = new TH1F("fTrackEffEta","Tracking efficienct vs eta",nbin,-1,1);
+  fFakeTrackEffEta = new TH1F("fFakeTrackEffEta","Efficiency for fake tracks vs eta",nbin,-1,1);
+
+  
+  fNGoodTracksPhi = new TH1F("fNGoodTracksPhi","Good tracks vs phi",nbin,0,6.28);
+  fNFoundTracksPhi = new TH1F("fNFoundTracksPhi","Found tracks vs phi",nbin,0,6.28);
+  fNFakeTracksPhi = new TH1F("fNFakeTracksPhi","Fake tracks vs phi",nbin,0,6.28);
+  fTrackEffPhi = new TH1F("fTrackEffPhi","Tracking efficienct vs phi",nbin,0,6.28);
+  fFakeTrackEffPhi = new TH1F("fFakeTrackEffPhi","Efficiency for fake tracks vs phi",nbin,0,6.28);
+  fFakesPt = new TProfile("fFakesPt","Number of ghosts vs pt",nbin,xlow,xup,0,100);
+  fFakesPhi = new TProfile("fFakesPhi","Number of ghosts vs phi",nbin,0,6.28,0,100);
+  fSlicesPt = new TProfile("fSlicesPt","Number of slices vs pt",nbin,xlow,xup,0,100);
+  fSlicesPhi = new TProfile("fSlicesPhi","Number of slices vs phi",nbin,0,6.28,0,100);
+  fPtRes2 = new TH1F("fPtRes2","Relative Pt resolution (Pt>2 GeV)",20,-100.,100.);
+  fEtaRes = new TH1F("fEtaRes","Eta resolution",30,-0.1,0.1);
+  fPsiRes = new TH1F("fPsiRes","Psi resolution",30,-0.1,0.1);
+
+  fNGoodTracksRad = new TH1F("fNGoodTracksRad","Good tracks vs distance(x,y)",1,0,0.2);    
+  fNFoundTracksRad = new TH1F("fNFoundTracksRad","Found tracks vs distance(x,y)",1,0,0.2);
+  fRadEff = new TH1F("fRadEff","Tracking efficiency vs distance(x,y)",1,0,0.2);
+  fNGoodTracksZ = new TH1F("fNGoodTracksZ","Good tracks vs distance(z)",1,0,0.2);    
+  fNFoundTracksZ = new TH1F("fNFoundTracksZ","Found tracks vs distance(z)",1,0,0.2);
+  fZEff = new TH1F("fZEff","Tracking efficiency vs distance(z)",1,0,0.2);
+}
+
+void Write2File()
+{
+  //Write histograms to file:
+  
+  TFile *of = TFile::Open("hough_efficiencies.root","RECREATE");
+  if(!of->IsOpen())
+    {
+      cout<<"Problems opening rootfile"<<endl;
+      return;
+    }
+  
+  of->cd();
+  fNtuppel->Write();
+  fPtRes->Write();
+  fNGoodTracksPt->Write();
+  fNFoundTracksPt->Write();
+  fNFakeTracksPt->Write();
+  fTrackEffPt->Write();
+  fFakeTrackEffPt->Write();
+  fNGoodTracksEta->Write();
+  fNFoundTracksEta->Write();
+  fNFakeTracksEta->Write();
+  fTrackEffEta->Write();
+  fFakeTrackEffEta->Write();
+
+  fNGoodTracksPhi->Write();
+  fNFoundTracksPhi->Write();
+  fNFakeTracksPhi->Write();
+  fTrackEffPhi->Write();
+  fFakeTrackEffPhi->Write();
+  
+  fFakesPt->Write();
+  fFakesPhi->Write();
+  fSlicesPt->Write();
+  fSlicesPhi->Write();
+
+  fPtRes2->Write();
+  fEtaRes->Write();
+  fPsiRes->Write();
+
+  fNGoodTracksRad->Write();
+  fNFoundTracksRad->Write(); 
+  fRadEff->Write();
+  fNGoodTracksZ->Write();
+  fNFoundTracksZ->Write();
+  fZEff->Write();
+
+  of->Close();
+}
+
+void CalcEffHistos()
+{  
+  fNFoundTracksPt->Sumw2(); fNGoodTracksPt->Sumw2();
+  fTrackEffPt->Divide(fNFoundTracksPt,fNGoodTracksPt,1,1.,"b");
+  fFakeTrackEffPt->Divide(fNFakeTracksPt,fNGoodTracksPt,1,1.,"b");
+  fTrackEffPt->SetMaximum(1.4);
+  fTrackEffPt->SetXTitle("P_{T} [GeV]");
+  fTrackEffPt->SetLineWidth(2);
+  fFakeTrackEffPt->SetFillStyle(3013);
+  fTrackEffPt->SetLineColor(4);
+  fFakeTrackEffPt->SetFillColor(2);
+
+  fNFoundTracksEta->Sumw2(); fNGoodTracksEta->Sumw2();
+  fTrackEffEta->Divide(fNFoundTracksEta,fNGoodTracksEta,1,1.,"b");
+  fFakeTrackEffEta->Divide(fNFakeTracksEta,fNGoodTracksEta,1,1.,"b");
+  fTrackEffEta->SetMaximum(1.4);
+  fTrackEffEta->SetXTitle("#eta");
+  fTrackEffEta->SetLineWidth(2);
+  fFakeTrackEffEta->SetFillStyle(3013);
+  fTrackEffEta->SetLineColor(4);
+  fFakeTrackEffEta->SetFillColor(2);
+       
+  fNFoundTracksPhi->Sumw2(); fNGoodTracksPhi->Sumw2();
+  fTrackEffPhi->Divide(fNFoundTracksPhi,fNGoodTracksPhi,1,1.,"b");
+  fFakeTrackEffPhi->Divide(fNFakeTracksPhi,fNGoodTracksPhi,1,1.,"b");
+  fTrackEffPhi->SetMaximum(1.4);
+  fTrackEffPhi->SetXTitle("#psi [rad]");
+  fTrackEffPhi->SetLineWidth(2);
+  fFakeTrackEffPhi->SetFillStyle(3013);
+  fTrackEffPhi->SetLineColor(4);
+  fFakeTrackEffPhi->SetFillColor(2);
+
+  fEtaRes->SetXTitle("#eta");
+  fPsiRes->SetXTitle("#psi [rad]");
+
+  fRadEff->Divide(fNFoundTracksRad,fNGoodTracksRad,1,1.,"b");
+  fRadEff->SetXTitle("Distance [cm]");
+  fZEff->Divide(fNFoundTracksZ,fNGoodTracksZ,1,1.,"b");
+  fZEff->SetXTitle("Distance [cm]");
+}
+
+void plotptres(Float_t ptmin)
+{
+  //Plot the relative pt resolution vs pt.
+  
+  const Int_t n=8;
+  
+  Double_t hltx[9] = {0.2,0.4,0.6,0.8,1.0,1.2,1.4,1.6,2.0};
+  Double_t hltxerr[n];
+  Double_t hltavx[n];
+  Double_t hlty[n];
+  Double_t hltyerr[n];
+  Char_t string[1024];
+
+  for(Int_t i=0; i<n; i++)
+    {
+      TH1F *hist = new TH1F("hist","",100,-50.*(float)(i+1),50.*(float)(i+1));
+      sprintf(string,"pt_gen > %f && pt_gen <= %f",hltx[i],hltx[i+1]);
+      fNtuppel->Draw("(pt_found-pt_gen)/pt_gen*100>>hist",string,"goff");
+      TF1 *f1 = new TF1("f1","gaus");
+      hist->Fit("f1","E");
+      hlty[i] = f1->GetParameter(2);
+      hltyerr[i] = f1->GetParError(2);
+      //hlty[i] = hist->GetRMS();
+      //hltyerr[i] = 0;
+      hltavx[i] = (hltx[i]+hltx[i+1])/2.0;
+      hltxerr[i] = (hltx[i+1]-hltx[i])/2.0;
+      cout<<string<<" "<<i<<" "<<hlty[i]<<" "<<hltyerr[i]<<endl;
+      delete f1;
+      delete hist;
+    }
+
+  sprintf(string,"pt_gen > %f",2.0);
+  fNtuppel->Draw("(pt_found-pt_gen)/pt_gen*100>>fPtRes2",string,"goff");
+  TF1 *f1 = new TF1("f1","gaus");
+  fPtRes2->Fit("f1","E");
+  fPtRes2->SetXTitle("#Delta P_{t} / P_{t} [%]");
+  fPtRes2->SetYTitle("#Delta N");
+  delete f1;
+
+  TGraphErrors *gr1 = new TGraphErrors(n,hltavx,hlty,hltxerr,hltyerr);
+  TCanvas *c1 = new TCanvas("c1","",1);
+  fPtRes = c1->DrawFrame(0.1,0,2,15);
+  gr1->SetTitle("");
+  gr1->Draw("P");
+  gr1->SetMarkerStyle(20);
+  gr1->SetLineWidth(2);
+  gr1->SetMarkerSize(1.3);
+  f1 = new TF1("f1","pol1",ptmin,2.0);
+  gr1->Fit("f1","ER");
+  fPtRes->SetXTitle("p_{t} [GeV]");
+  fPtRes->SetYTitle("#Delta P_{t} / P_{t} [%]");
+  delete f1;
+}
index e8f80fc..7e84cf3 100644 (file)
    and written to file.
 */
 
-void evaltracker(char *path,char *trackpath,char *offlinepath,int nevent=1)
+#ifndef __CINT__
+#include "AliL3Logger.h"
+#include "AliL3Evaluate.h"
+#include "AliL3FileHandler.h"
+#include "AliL3DigitData.h"
+#include "AliL3Transform.h"
+#include "AliLevel3.h"
+#include <TROOT.h>
+#include <TNtuple.h>
+#include <TRandom.h>
+#include <TSystem.h>
+#include <TStyle.h>
+#include <TFile.h>
+#include <TGraphErrors.h>
+#include <TCanvas.h>
+#include <TF1.h>
+#include <stdio.h>
+#include <iostream.h>
+#include <time.h>
+#endif
+
+void evaltracker(Char_t *path="./",Char_t *trackpath="./tracker/",char *offlinepath="./",int nevent=1)
 {
   //Make sure you got the correct parameters:
   AliL3Transform::Init(path,kTRUE);
   
   //Define which slices to include:
-  int slicerange[2]={0,35};
+  Int_t slicerange[2]={0,35};
   
   //Define which padrows to include (should normally be all)
-  int rowrange[2]={AliL3Transform::GetFirstRow(-1),AliL3Transform::GetLastRow(-1)};
+  //Int_t rowrange[2]={AliL3Transform::GetFirstRow(-1),AliL3Transform::GetLastRow(-1)};
   
   //Define the minimum number of clusters on a simulated track to be found:
-  int good = (int)(0.4*AliL3Transform::GetNRows());
+  Int_t good = (Int_t)(0.4*AliL3Transform::GetNRows());
   
   //Define the minumum number of clusters on a found track to be found (should normally be the same as above)
-  int nclusters = (int)(0.4*AliL3Transform::GetNRows());
+  Int_t nclusters = (Int_t)(0.4*AliL3Transform::GetNRows());
   
   //Define which pt range to include
-  float ptmin = 0.1;
-  float ptmax = 4.;
+  Float_t ptmin = 0.1;
+  Float_t ptmax = 4.;
   
   //Define the maximum ratio of false clusters which are allowed on a found track (default=0.1 -> 10%)
-  float maxfalseratio = 0.1;
-  
-  a = new AliL3Evaluate(trackpath,nclusters,good,ptmin,ptmax,slicerange);
+  Float_t maxfalseratio = 0.1;
   
+  AliL3Evaluate *a = new AliL3Evaluate(trackpath,nclusters,good,ptmin,ptmax,slicerange);
   a->CreateHistos(20,0.1,4.1);//(nbins in pt,minpt,maxpt) -> the same as used by standard offline
   a->SetMaxFalseClusters(maxfalseratio);
   //a->SetStandardComparison(kFALSE); //use AliTPCComparison_HLT
 
   //Loop over the number of events:
-  for(int event=0; event<nevent; event++)
+  for(Int_t event=0; event<nevent; event++)
     {
       cout<<"Processing event "<<event<<endl;      
       
@@ -73,24 +93,20 @@ void evaltracker(char *path,char *trackpath,char *offlinepath,int nevent=1)
   
 }
 
-void ploteff(char *rootfile="hlt_efficiencies.root")
+void ploteff(Char_t *rootfile="hlt_efficiencies.root")
 {
   //Plot the efficiency vs pt.
   
   gStyle->SetOptStat(0);
-  char filename[1024];
-  
-  double hltx[22];
-  double hlty[22];
-  double hltxerr[22];
-  double hltyerr[22];
-  int i;
   
-  file = TFile::Open(rootfile);
-  
-  int hltcounter=0;
-  h = (TH1*)fTrackEffPt;
-  for(i=0; i<22; i++)
+  Double_t hltx[22];
+  Double_t hlty[22];
+  Double_t hltxerr[22];
+  Double_t hltyerr[22];
+  TFile *file = TFile::Open(rootfile);
+  Int_t hltcounter=0;
+  TH1* h = (TH1*)file->Get("fTrackEffPt");
+  for(Int_t i=0; i<22; i++)
     {
       if(h->GetBinContent(i)==0) continue;
       hltx[hltcounter] = h->GetBinCenter(i);
@@ -101,10 +117,10 @@ void ploteff(char *rootfile="hlt_efficiencies.root")
     }
   file->Close();
   
-  gr1 = new TGraphErrors(hltcounter,hltx,hlty,hltxerr,hltyerr);
+  TGraphErrors *gr1 = new TGraphErrors(hltcounter,hltx,hlty,hltxerr,hltyerr);
   
-  c1 = new TCanvas("c1","",1);
-  hist = c1->DrawFrame(0.1,0,3,1.4);
+  TCanvas *c1 = new TCanvas("c1","",1);
+  TH1F *hist = c1->DrawFrame(0.1,0,3,1.4);
   c1->SetGridx();
   c1->SetGridy();
   gr1->SetTitle("");
@@ -114,31 +130,28 @@ void ploteff(char *rootfile="hlt_efficiencies.root")
   gr1->SetMarkerStyle(20);
   gr1->SetLineWidth(2);
   gr1->SetMarkerSize(1.3);
-  
 }
 
-void plotptres(char *rootfile="hlt_efficiencies.root")
+void plotptres(Char_t *rootfile="hlt_efficiencies.root")
 {
   //Plot the relative pt resolution vs pt.
   
-  const int n=5;
+  const Int_t n=15;
   
-  double hltx[15] = {0.2,0.4,0.6,0.8,1.0,1.2,1.4,1.6,1.8,2.0,2.2,2.4,2.6,2.8,3.0};
-  double hlty[n];
-  double hltyerr[n];
-  int i;
+  Double_t hltx[15] = {0.2,0.4,0.6,0.8,1.0,1.2,1.4,1.6,1.8,2.0,2.2,2.4,2.6,2.8,3.0};
+  Double_t hlty[n];
+  Double_t hltyerr[n];
+  Char_t string[1024];
   
-  char string[1024];
+  TFile *file = TFile::Open(rootfile);
+  TH1F *hist = new TH1F("hist","",100,-10,10);
+  TNtuple *fNtuppel=(TNtuple*)file->Get("fNtuppel");
   
-  file = TFile::Open(rootfile);
-  
-  hist = new TH1F("hist","",100,-10,10);
-  
-  for(i=0; i<n; i++)
+  for(Int_t i=0; i<n; i++)
     {
       sprintf(string,"pt_gen > %f && pt_gen <= %f && nHits>63",hltx[i]-0.1,hltx[i]+0.1);
       fNtuppel->Draw("(pt_found-pt_gen)/pt_gen*100>>hist",string,"goff");
-      f1 = new TF1("f1","gaus");
+      TF1 *f1 = new TF1("f1","gaus");
       hist->Fit("f1","QN");
       hlty[i] = f1->GetParameter(2);
       hltyerr[i] = f1->GetParError(2);
@@ -148,10 +161,9 @@ void plotptres(char *rootfile="hlt_efficiencies.root")
   delete hist;
   file->Close();
 
-  gr1 = new TGraphErrors(n,hltx,hlty,0,hltyerr);
-  
-  c1 = new TCanvas("c1","",1);
-  hist = c1->DrawFrame(0.1,0,3,5);
+  TGraphErrors *gr1 = new TGraphErrors(n,hltx,hlty,0,hltyerr);
+  TCanvas *c1 = new TCanvas("c1","",1);
+  hist = c1->DrawFrame(0.1,0,3,15);
   c1->SetGridx();
   c1->SetGridy();
   gr1->SetTitle("");
@@ -161,5 +173,4 @@ void plotptres(char *rootfile="hlt_efficiencies.root")
   gr1->SetMarkerSize(1.3);
   hist->SetXTitle("p_{t} [GeV]");
   hist->SetYTitle("#Delta P_{T} / P_{T} [%]");
-  
 }
index 859402f..eef4c34 100644 (file)
@@ -62,6 +62,7 @@
       gSystem->Load("MakePileup_C.so");
       gSystem->Load("Read_C.so");
       gSystem->Load("runhough_C.so");
+      gSystem->Load("runrowhough_C.so");
       gSystem->Load("deconvclusters_C.so");
       gSystem->Load("runtracker_C.so");
       gErrorIgnoreLevel=saveErrIgLevel;
diff --git a/HLT/exa/runrowhough.C b/HLT/exa/runrowhough.C
new file mode 100644 (file)
index 0000000..22996d0
--- /dev/null
@@ -0,0 +1,159 @@
+// $Id$
+
+#ifndef __CINT__
+#include "AliL3Logger.h"
+#include "AliL3FileHandler.h"
+#include "AliL3DigitData.h"
+#include "AliL3Transform.h"
+#include "AliL3Hough.h"
+#include "AliL3TrackArray.h"
+#include "AliL3Track.h"
+#include "AliL3HoughTrack.h"
+#include "AliL3Fitter.h"
+#include "AliL3ClusterFitter.h"
+#include "AliL3Vertex.h"
+#include "AliL3Benchmark.h"
+#include <AliRunLoader.h>
+#include <AliStack.h>
+#include <TParticle.h>
+#include <TNtuple.h>
+#include <TRandom.h>
+#include <TSystem.h>
+#include <TStopwatch.h>
+#include <TBenchmark.h>
+#include <stdio.h>
+#include <iostream.h>
+#include <time.h>
+#endif
+
+Int_t runrowhough(Char_t *path="./",Char_t *outpath="./fitter",int s1=0,int s2=35,int nevent=1,Bool_t skip=kTRUE)
+{
+  Bool_t isinit=AliL3Transform::Init(path,kTRUE);
+  if(!isinit){
+    cerr << "Could not create transform settings, please check log for error messages!" << endl;
+    return 1;
+  }
+  Float_t ptmin = 0.1*AliL3Transform::GetSolenoidField();
+  Float_t zvertex;
+
+  {
+    AliRunLoader *rl = AliRunLoader::Open("galice.root");
+    rl->LoadHeader();
+    rl->LoadKinematics();
+    AliStack* stack = rl->Stack();
+    TParticle *orig = (TParticle*)stack->Particle(0);
+    Float_t xori = orig->Vx();
+    Float_t yori = orig->Vy();
+    zvertex = orig->Vz();
+    cout<<" Primary vertex at ("<<xori<<","<<yori<<","<<zvertex<<")"<<endl; 
+    if (rl->LoadgAlice()) {
+      cerr<<"Error occured while loading gAlice"<<endl;
+      return 1;
+    }
+    delete rl;
+  }
+
+  cout<<" Hough Tranform will run with ptmin="<<ptmin<<" and zvertex="<<zvertex<<endl;
+  
+  AliL3Benchmark *fBenchmark = new AliL3Benchmark();
+  AliL3Hough *hough = new AliL3Hough();
+  hough->SetThreshold(4);
+  hough->SetTransformerParams(140,76,ptmin,-1);
+  hough->SetPeakThreshold(50,-1);
+  hough->Init(path, kFALSE, 100, kFALSE,4,0,0,zvertex);
+  hough->SetAddHistograms();
+
+  for(int ev=0; ev<nevent; ev++)
+    {
+      for(int slice=s1; slice<=s2; slice++)
+       {
+         cout<<"Processing slice "<<slice<<endl;
+         hough->ReadData(slice,ev);
+         hough->Transform();
+         hough->AddAllHistogramsRows();
+         hough->FindTrackCandidates();
+         hough->AddTracks();
+       }
+      hough->WriteTracks(outpath);
+      Char_t bname[100];
+      sprintf(bname,"rowhough_%d",ev);
+      hough->DoBench(bname);
+
+      if(!skip) {
+       // Run cluster fitter
+       AliL3ClusterFitter *fitter = new AliL3ClusterFitter(path);
+
+       // Set debug flag for the cluster fitter
+       //  fitter->Debug();
+
+       // Setting fitter parameters
+       fitter->SetInnerWidthFactor(1,1.5);
+       fitter->SetOuterWidthFactor(1,1.5);
+       fitter->SetNmaxOverlaps(5);
+  
+       //  fitter->SetChiSqMax(5,kFALSE); //isolated clusters
+       fitter->SetChiSqMax(5,kTRUE);  //overlapping clusters
+
+       Int_t rowrange[2] = {0,AliL3Transform::GetNRows()-1};
+
+       // Takes input from global hough tracks produced by HT
+       fitter->LoadSeeds(rowrange,kFALSE,ev,zvertex);
+
+       UInt_t ndigits;
+
+       for(int slice=s1; slice<=s2; slice++)
+         {
+           for(Int_t ipatch = 0; ipatch < AliL3Transform::GetNPatches(); ipatch++)
+             {
+               // Read digits
+               hough->GetMemHandler(ipatch)->Free();
+               hough->GetMemHandler(ipatch)->Init(slice,ipatch);
+               AliL3DigitRowData *digits = (AliL3DigitRowData *)hough->GetMemHandler(ipatch)->AliAltroDigits2Memory(ndigits,ev);
+
+               fBenchmark->Start("Fitter Init");
+               fitter->Init(slice,ipatch);
+               fBenchmark->Stop("Fitter Init");
+               fitter->SetInputData(digits);
+               fBenchmark->Start("Fitter cluster finder");
+               fitter->FindClusters();
+               fBenchmark->Stop("Fitter cluster finder");
+               fitter->WriteClusters();
+             }
+         }
+
+       // Refit of the clusters
+       AliL3Vertex vertex;
+       //The seeds are the input tracks from circle HT
+       AliL3TrackArray *tracks = fitter->GetSeeds();
+       AliL3Fitter *ft = new AliL3Fitter(&vertex,1);
+
+       ft->LoadClusters("./fitter/",ev,kFALSE);
+       fBenchmark->Start("Track fitter");
+       for(Int_t i=0; i<tracks->GetNTracks(); i++)
+         {
+           AliL3Track *track = tracks->GetCheckedTrack(i);
+           if(!track) continue;
+           if(track->GetNHits() < 20) continue;
+           ft->SortTrackClusters(track);
+           ft->FitHelix(track);
+           track->UpdateToFirstPoint();
+         }
+       fBenchmark->Stop("Track fitter");
+       delete ft;
+        
+       //Write the final tracks
+       fitter->WriteTracks(20);
+
+       delete fitter;
+      }
+    }
+
+  if(!skip) {
+    fBenchmark->Analyze("fitter");
+  }
+
+  hough->DoBench("rowhough");
+  delete hough;
+  delete fBenchmark;
+  return 0;
+}
index a10f6bd..ffd5af5 100644 (file)
@@ -103,8 +103,9 @@ void runtracker(Int_t minslice=0,Int_t maxslice=35,Char_t* path="./",Int_t neven
       //a->DoMc();     /*do monte carlo identification*/
       a->WriteFiles(opath); /*enable output*/
       a->ProcessEvent(minslice,maxslice);
-      a->DoBench("benchmark_0");
-
+      Char_t bname[100];
+      sprintf(bname,"benchmark_tracker_%d",ev);
+      a->DoBench(bname);
       delete a;
     } // event loop
 }
index 3c2ef5f..ccc9cee 100644 (file)
@@ -8,14 +8,22 @@
 #include "AliL3Logging.h"
 #include "AliL3Histogram.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
+/** \class AliL3Histogram
+<pre>
 //_____________________________________________________________
 // AliL3Histogram
 //
 // 2D histogram class
+//
+</pre>
+*/
+
+//uncomment if you want overflow checks
+//#define _IFON_
 
 ClassImp(AliL3Histogram)
 
@@ -296,6 +304,30 @@ Double_t AliL3Histogram::GetBinCenterY(Int_t ybin) const
   return fYmin + (ybin-0.5) * fBinwidthY;
 }
 
+Double_t AliL3Histogram::GetPreciseBinCenterX(Float_t xbin) const
+{
+  if(xbin < fFirstXbin || xbin > fLastXbin)
+    {
+      LOG(AliL3Log::kError,"AliL3Histogram::GetBinCenterX","xbin")
+       <<"Bin-value out of range "<<xbin<<ENDLOG;
+      return -1;
+    }
+  //  return fXmin + (xbin-1) * fBinwidthX + 0.5*fBinwidthX;
+  return fXmin + (xbin-0.5) * fBinwidthX;
+}
+
+Double_t AliL3Histogram::GetPreciseBinCenterY(Float_t ybin) const
+{
+  if(ybin < fFirstYbin || ybin > fLastYbin)
+    {
+      LOG(AliL3Log::kError,"AliL3Histogram::GetBinCenterY","ybin")
+       <<"Bin-value out of range "<<ybin<<ENDLOG;
+      return -1;
+    }
+  //  return fYmin + (ybin-1) * fBinwidthY + 0.5*fBinwidthY;
+  return fYmin + (ybin-0.5) * fBinwidthY;
+}
+
 void AliL3Histogram::Draw(Char_t *option)
 {
 #ifdef use_root
index 339fc40..79b4aa0 100644 (file)
@@ -77,6 +77,8 @@ class AliL3Histogram {
   Double_t GetYmax() const {return fYmax;}
   virtual Double_t GetBinCenterX(Int_t xbin) const;
   virtual Double_t GetBinCenterY(Int_t ybin) const;
+  Double_t GetPreciseBinCenterX(Float_t xbin) const;
+  Double_t GetPreciseBinCenterY(Float_t ybin) const;
   Double_t GetBinWidthX() const {return fBinwidthX;}
   Double_t GetBinWidthY() const {return fBinwidthY;}
   Int_t GetFirstXbin() const {return fFirstXbin;}
index add80d7..2bc348e 100644 (file)
@@ -8,7 +8,7 @@
 #include "AliL3Logging.h"
 #include "AliL3Histogram1D.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 8bb7c85..e9bf95d 100644 (file)
@@ -9,7 +9,7 @@
 #include "AliL3Transform.h"
 #include "AliL3Track.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 57b27b1..5cb48ba 100644 (file)
@@ -16,7 +16,7 @@
 #include "AliL3HoughClusterTransformer.h"
 #include "AliL3HoughTransformerLUT.h"
 #include "AliL3HoughTransformerVhdl.h"
-#include "AliL3HoughTransformerGap.h"
+#include "AliL3HoughTransformerRow.h"
 #include "AliL3HoughMaxFinder.h"
 #include "AliL3Benchmark.h"
 #ifdef use_aliroot
@@ -32,7 +32,7 @@
 #include "AliL3HoughTrack.h"
 #include "AliL3DDLDataFileHandler.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
@@ -51,6 +51,7 @@ using namespace std;
 // hough->FindTrackCandidates();
 // 
 // AliL3TrackArray *tracks = hough->GetTracks(patch);
+//
 //</pre>
 */
 
@@ -83,8 +84,10 @@ AliL3Hough::AliL3Hough()
   fCurrentSlice     = 0;
   fEvent            = 0;
   
-  fKappaSpread=6;
-  fPeakRatio=0.5;
+  fKappaSpread = 6;
+  fPeakRatio   = 0.5;
+  fInputFile   = 0;
+  fInputPtr    = 0;
   
   SetTransformerParams();
   SetThreshold();
@@ -96,7 +99,7 @@ AliL3Hough::AliL3Hough()
 #endif
 }
 
-AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv,Char_t *infile)
+AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr)
 {
   fBinary = binary;
   strcpy(fPath,path);
@@ -108,11 +111,20 @@ AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bi
   fVersion       = tv;
   fKappaSpread=6;
   fPeakRatio=0.5;
-  if(!fBinary)
-    fInputFile = infile;
-  else
+  if(!fBinary) {
+    if(infile) {
+      fInputFile = infile;
+      fInputPtr = 0;
+    }
+    else {
+      fInputFile = 0;
+      fInputPtr = ptr;
+    }
+  }
+  else {
     fInputFile = 0;
-
+    fInputPtr = 0;
+  }
 #ifdef use_aliroot
   //just be sure that index is empty for new event
     AliL3FileHandler::CleanStaticIndex(); 
@@ -170,7 +182,7 @@ void AliL3Hough::CleanUp()
   //cout << "Cleaned class mem " << endl;
 }
 
-void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv,Char_t *infile)
+void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr,Float_t zvertex)
 {
   fBinary = binary;
   strcpy(fPath,path);
@@ -178,10 +190,21 @@ void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit
   fWriteDigits  = kFALSE;
   fUse8bits     = bit8;
   fVersion      = tv;
-  if(!fBinary)
-    fInputFile = infile;
-  else
+  if(!fBinary) {
+    if(infile) {
+      fInputFile = infile;
+      fInputPtr = 0;
+    }
+    else {
+      fInputFile = 0;
+      fInputPtr = ptr;
+    }
+  }
+  else {
     fInputFile = 0;
+    fInputPtr = 0;
+  }
+  fZVertex = zvertex;
 
   Init(); //do the rest
 }
@@ -214,7 +237,7 @@ void AliL3Hough::Init(Bool_t doit, Bool_t addhists)
        fHoughTransformer[i] = new AliL3HoughTransformerVhdl(0,i,fNEtaSegments,fNSaveIterations);
        break;
       case 4:
-       fHoughTransformer[i] = new AliL3HoughTransformerGap(0,i,fNEtaSegments);
+       fHoughTransformer[i] = new AliL3HoughTransformerRow(0,i,fNEtaSegments,kFALSE,fZVertex);
        break;
       default:
        fHoughTransformer[i] = new AliL3HoughTransformer(0,i,fNEtaSegments,kFALSE,kFALSE);
@@ -237,12 +260,19 @@ void AliL3Hough::Init(Bool_t doit, Bool_t addhists)
 #ifdef use_aliroot
        {
          if(!fInputFile) {
-           /* In case of reading digits file */
-           fMemHandler[i] = new AliL3FileHandler(kTRUE); //use static index
-           if(!fBinary) {
-             Char_t filename[1024];
-             sprintf(filename,"%s/digitfile.root",fPath);
-              fMemHandler[i]->SetAliInput(filename);
+           if(!fInputPtr) {
+             /* In case of reading digits file */
+             fMemHandler[i] = new AliL3FileHandler(kTRUE); //use static index
+             if(!fBinary) {
+               Char_t filename[1024];
+               sprintf(filename,"%s/digitfile.root",fPath);
+               fMemHandler[i]->SetAliInput(filename);
+             }
+           }
+           else {
+             /* In case of reading from DATE */
+             fMemHandler[i] = new AliL3DDLDataFileHandler();
+             fMemHandler[i]->SetReaderInput(fInputPtr,-1);
            }
          }
          else {
@@ -374,8 +404,12 @@ void AliL3Hough::Process(Int_t minslice,Int_t maxslice)
     {
       ReadData(i);
       Transform();
-      if(fAddHistograms)
-       AddAllHistograms();
+      if(fAddHistograms) {
+       if(fVersion != 4)
+         AddAllHistograms();
+       else
+         AddAllHistogramsRows();
+      }
       FindTrackCandidates();
       //Evaluate();
       //fGlobalMerger->FillTracks(fTracks[0],i);
@@ -438,7 +472,9 @@ void AliL3Hough::Transform(Int_t *row_range)
   initTime = GetCpuTime();
   for(Int_t i=0; i<fNPatches; i++)
     {
-      fHoughTransformer[i]->Reset();//Reset the histograms
+      // In case of Row transformer reset the arrays only once
+      if((fVersion != 4) || (i == 0)) 
+       fHoughTransformer[i]->Reset();//Reset the histograms
       fBenchmark->Start("Hough Transform");
       if(!row_range)
        fHoughTransformer[i]->TransformCircle();
@@ -590,6 +626,45 @@ void AliL3Hough::AddAllHistograms()
     <<"Adding histograms in "<<cpuTime*1000<<" ms"<<ENDLOG;
 }
 
+void AliL3Hough::AddAllHistogramsRows()
+{
+  //Add the histograms within one etaslice.
+  //Resulting histogram are in patch=0.
+
+  Double_t initTime,cpuTime;
+  initTime = GetCpuTime();
+  fBenchmark->Start("Add HistogramsRows");
+
+  for(Int_t i=0; i<fNEtaSegments; i++)
+    {
+      UChar_t *rowcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetRowCount(i);
+      UChar_t *gapcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetGapCount(i);
+
+      AliL3Histogram *hist = fHoughTransformer[0]->GetHistogram(i);
+      Int_t xmin = hist->GetFirstXbin();
+      Int_t xmax = hist->GetLastXbin();
+      Int_t ymin = hist->GetFirstYbin();
+      Int_t ymax = hist->GetLastYbin();
+
+      for(Int_t ybin=ymin; ybin<=ymax; ybin++)
+       {
+         for(Int_t xbin=xmin; xbin<=xmax; xbin++)
+           {
+             Int_t bin = hist->GetBin(xbin,ybin);
+             if(gapcount[bin] < 3)
+               //              hist->AddBinContent(bin,rowcount[bin]);
+               hist->AddBinContent(bin,rowcount[bin]/gapcount[bin]);
+           }
+       }
+    }
+
+  fBenchmark->Stop("Add HistogramsRows");
+  fAddHistograms = kTRUE;
+  cpuTime = GetCpuTime() - initTime;
+  LOG(AliL3Log::kInformational,"AliL3Hough::AddAllHistogramsRows()","Timing")
+    <<"Adding histograms in "<<cpuTime*1000<<" ms"<<ENDLOG;
+}
+
 void AliL3Hough::AddTracks()
 {
   if(!fTracks[0])
@@ -636,9 +711,12 @@ void AliL3Hough::FindTrackCandidates()
          fPeakFinder->SetHistogram(hist);
 
          fPeakFinder->SetThreshold(fPeakThreshold[i]);
-         fPeakFinder->FindAdaptedPeaks(fKappaSpread,fPeakRatio);
-
-         //fPeakFinder->FindMaxima(fPeakThreshold[i]); //Simple maxima finder
+         if(fVersion != 4) {
+           fPeakFinder->FindAdaptedPeaks(fKappaSpread,fPeakRatio);
+           //fPeakFinder->FindMaxima(fPeakThreshold[i]); //Simple maxima finder
+         }
+         else /* row transformer needs different peak finder method */
+           fPeakFinder->FindAdaptedRowPeaks(1,2,0);
          
          for(Int_t k=0; k<fPeakFinder->GetEntries(); k++)
            {
@@ -647,6 +725,11 @@ void AliL3Hough::FindTrackCandidates()
              track->SetEtaIndex(j);
              track->SetEta(tr->GetEta(j,fCurrentSlice));
              track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5));
+#ifdef do_mc
+             Int_t label = tr->GetTrackID(j,fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k));
+             track->SetMCid(label);
+             cout<<"Track found with label "<<label<<" at "<<fPeakFinder->GetXPeak(k)<<" "<<fPeakFinder->GetYPeak(k)<<" with weight "<<fPeakFinder->GetWeight(k)<<endl; 
+#endif
            }
        }
       cout<<"Found "<<fTracks[i]->GetNTracks()<<" tracks in patch "<<i<<endl;
index 88e7de7..519e3a0 100644 (file)
@@ -21,7 +21,7 @@ class AliL3Hough {
   
  private:
   Char_t *fInputFile;//!
-
+  Char_t *fInputPtr;//!
   Char_t fPath[1024];
   Bool_t fBinary;
   Bool_t fAddHistograms;
@@ -48,6 +48,8 @@ class AliL3Hough {
   Int_t fKappaSpread;
   Float_t fPeakRatio;
 
+  Float_t fZVertex;
+
   AliL3MemHandler **fMemHandler; //!
   AliL3HoughBaseTransformer **fHoughTransformer; //!
   AliL3HoughEval **fEval; //!
@@ -65,10 +67,10 @@ class AliL3Hough {
  public:
   
   AliL3Hough(); 
-  AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments=100,Bool_t bit8=kFALSE,Int_t tv=0,Char_t *infile=0);
+  AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments=100,Bool_t bit8=kFALSE,Int_t tv=0,Char_t *infile=0,Char_t *ptr=0);
   virtual ~AliL3Hough();
   
-  void Init(Char_t *path,Bool_t binary,Int_t n_eta_segments=100,Bool_t bit8=kFALSE,Int_t tv=0,Char_t *infile=0);
+  void Init(Char_t *path,Bool_t binary,Int_t n_eta_segments=100,Bool_t bit8=kFALSE,Int_t tv=0,Char_t *infile=0,Char_t *ptr=0,Float_t zvertex=0.0);
   void Init(Bool_t doit=kFALSE, Bool_t addhists=kFALSE);
 
   void Process(Int_t minslice,Int_t maxslice);
@@ -82,6 +84,7 @@ class AliL3Hough {
 
   void FindTrackCandidates();
   void AddAllHistograms();
+  void AddAllHistogramsRows();
   Int_t Evaluate(Int_t road_width=1,Int_t nrowstomiss=1);
   void EvaluatePatch(Int_t i,Int_t road_width,Int_t nrowstomiss);
   void WriteTracks(Int_t slice,Char_t *path="./");
index 74d2efe..1e11850 100644 (file)
@@ -11,7 +11,8 @@
 #include "AliL3Histogram.h"
 #include "AliL3HoughBaseTransformer.h"
 
-
+/** \class AliL3HoughBaseTransformer
+<pre>
 //_____________________________________________________________
 // AliL3HoughBaseTransformer
 //
@@ -19,6 +20,9 @@
 //
 // This is an abstract class, and is only meant to provide the interface
 // to the different implementations.
+//
+</pre>
+*/
 
 ClassImp(AliL3HoughBaseTransformer)
 
@@ -34,9 +38,10 @@ AliL3HoughBaseTransformer::AliL3HoughBaseTransformer()
   fEtaMax = 0;
   fLowerThreshold = 0;
   fUpperThreshold = 1023;
+  fZVertex = 0.0;
 }
 
-AliL3HoughBaseTransformer::AliL3HoughBaseTransformer(Int_t slice,Int_t patch,Int_t n_eta_segments)
+AliL3HoughBaseTransformer::AliL3HoughBaseTransformer(Int_t slice,Int_t patch,Int_t n_eta_segments,Float_t zvertex)
 {
   fDigitRowData = 0;
 
@@ -47,6 +52,7 @@ AliL3HoughBaseTransformer::AliL3HoughBaseTransformer(Int_t slice,Int_t patch,Int
   fEtaMax = 0;
   fLowerThreshold = 3;
   fUpperThreshold = 1023;
+  fZVertex = zvertex;
 
   Init(slice,patch,n_eta_segments);
 }
index 411ac7d..03b48b4 100644 (file)
@@ -6,10 +6,15 @@
 #include "AliL3RootTypes.h"
 
 #ifdef do_mc
-const UInt_t MaxTrack=35;
+const UInt_t MaxTrack=100;
 struct TrackIndex {
   Int_t fLabel[MaxTrack];
+#ifdef ROWHOUGH
+  UChar_t fNHits[MaxTrack];
+  UChar_t fCurrentRow[MaxTrack];
+#else
   Int_t fNHits[MaxTrack];
+#endif
 };
 typedef struct TrackIndex TrackIndex;
 #endif
@@ -30,11 +35,13 @@ class AliL3HoughBaseTransformer {
   Int_t fUpperThreshold;
   
   AliL3DigitRowData *fDigitRowData; //!
-  
+
+  Float_t fZVertex;
+
  public:
 
   AliL3HoughBaseTransformer(); 
-  AliL3HoughBaseTransformer(Int_t slice,Int_t patch,Int_t n_eta_segments);
+  AliL3HoughBaseTransformer(Int_t slice,Int_t patch,Int_t n_eta_segments,Float_t zvertex=0.0);
   virtual ~AliL3HoughBaseTransformer();
   
   void SetInputData(UInt_t ndigits,AliL3DigitRowData *ptr) {fDigitRowData = ptr;}
@@ -66,7 +73,8 @@ class AliL3HoughBaseTransformer {
   Int_t GetUpperThreshold() {return fUpperThreshold;}
   Double_t GetEtaMin() {return fEtaMin;}
   Double_t GetEtaMax() {return fEtaMax;}
-  
+  Float_t GetZVertex() {return fZVertex;}
+
   AliL3DigitRowData *GetDataPointer() {return fDigitRowData;}
  
   virtual Int_t GetEtaIndex(Double_t eta) = 0;
index 4ba4a55..7020dd6 100644 (file)
@@ -14,7 +14,7 @@
 #include "AliL3Histogram.h"
 #include "AliL3ClustFinderNew.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 40e42b5..6f3fe9f 100644 (file)
@@ -21,7 +21,7 @@
 #include "AliL3MemHandler.h"
 #include "AliL3HoughDisplay.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 536a061..48f207e 100644 (file)
 #include "AliL3Histogram.h"
 #include "AliL3Histogram1D.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
+/** /class AliL3HoughEval
+//<pre>
 //_____________________________________________________________
 // AliL3HoughEval
 //
 // Evaluation class for tracklets produced by the Hough transform.
+//
+</pre>
+*/
 
 ClassImp(AliL3HoughEval)
 
@@ -63,6 +68,7 @@ void AliL3HoughEval::InitTransformer(AliL3HoughBaseTransformer *transformer)
   fNEtaSegments = fHoughTransformer->GetNEtaSegments();
   fEtaMin = fHoughTransformer->GetEtaMin();
   fEtaMax = fHoughTransformer->GetEtaMax();
+  fZVertex = fHoughTransformer->GetZVertex();
   GenerateLUT();
 }
 
@@ -302,6 +308,7 @@ void AliL3HoughEval::DisplayEtaSlice(Int_t eta_index,AliL3Histogram *hist)
          Int_t sector,row;
          AliL3Transform::Slice2Sector(fSlice,padrow,sector,row);
          AliL3Transform::Raw2Local(xyz,sector,row,pad,time);
+         xyz[2] -= fZVertex;
          Double_t eta = AliL3Transform::GetEta(xyz);
          Int_t pixel_index = fHoughTransformer->GetEtaIndex(eta);//(Int_t)(eta/etaslice);
          if(pixel_index != eta_index) continue;
index 44a98ac..553f80c 100644 (file)
@@ -25,6 +25,7 @@ class AliL3HoughEval {
   Int_t fNumOfPadsToLook;
   Int_t fNumOfRowsToMiss;
   AliL3Histogram1D **fEtaHistos; //!
+  Float_t fZVertex;
 
   //Flags
   Bool_t fRemoveFoundTracks;
@@ -53,6 +54,7 @@ class AliL3HoughEval {
   void SetNumOfRowsToMiss(Int_t i) {fNumOfRowsToMiss = i;}
   void SetNumOfPadsToLook(Int_t i) {fNumOfPadsToLook = i;}
   void SetSlice(Int_t i) {fSlice=i;}
+  void SetZVertex(Float_t zvertex) {fZVertex=zvertex;}
 
   ClassDef(AliL3HoughEval,1) //Hough transform verfication class
 
index 7010277..96f75e8 100644 (file)
@@ -12,7 +12,7 @@
 #include "AliL3Transform.h"
 #include "AliL3TrackArray.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index f4a52db..1629cfb 100644 (file)
@@ -12,7 +12,7 @@
 #pragma link C++ class AliL3HoughTransformerLUT;
 #pragma link C++ class AliL3HoughTransformerVhdl;
 #pragma link C++ class AliL3HoughTransformerNew;
-#pragma link C++ class AliL3HoughTransformerGap;
+#pragma link C++ class AliL3HoughTransformerRow;
 #ifndef Darwin
 #pragma link C++ class AliL3HoughTrack;
 #endif
@@ -35,4 +35,3 @@
 #endif
 
 #endif
-
index 488f194..02104f9 100644 (file)
 #include "AliL3TrackArray.h"
 #include "AliL3HoughTrack.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
+/** \class AliL3HoughMaxFinder
+<pre>
 //_____________________________________________________________
 // AliL3HoughMaxFinder
 //
 // Maximum finder
+//
+</pre>
+*/
 
 ClassImp(AliL3HoughMaxFinder)
 
@@ -44,7 +49,6 @@ AliL3HoughMaxFinder::AliL3HoughMaxFinder()
 #endif
 }
 
-
 AliL3HoughMaxFinder::AliL3HoughMaxFinder(Char_t *histotype,Int_t nmax,AliL3Histogram *hist)
 {
   //Constructor
@@ -68,7 +72,6 @@ AliL3HoughMaxFinder::AliL3HoughMaxFinder(Char_t *histotype,Int_t nmax,AliL3Histo
   fThreshold=0;
 }
 
-
 AliL3HoughMaxFinder::~AliL3HoughMaxFinder()
 {
   //Destructor
@@ -211,7 +214,6 @@ void AliL3HoughMaxFinder::FindBigMaxima()
                  value[bin_index]=hist->GetBinContent(bin[bin_index]);
                  bin_index++;
                }
-             
            }
          if(value[12]==0) continue;
          Int_t b=0;
@@ -369,12 +371,13 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio)
   
   if(hist->GetNEntries() == 0)
     return;
-  
+
   Int_t xmin = hist->GetFirstXbin();
   Int_t xmax = hist->GetLastXbin();
   Int_t ymin = hist->GetFirstYbin();
   Int_t ymax = hist->GetLastYbin();
-  
+
+
   //Start by looking for pre-peaks:
   
   Window **local_maxima = new Window*[hist->GetNbinsY()];
@@ -618,6 +621,189 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio)
   delete [] maxs;
 }
 
+struct PreYPeak 
+{
+  Int_t start_position;
+  Int_t end_position;
+  Int_t value;
+  Int_t prev_value;
+};
+
+struct Pre2DPeak 
+{
+  Int_t start_x_position;
+  Int_t end_x_position;
+  Int_t start_y_position;
+  Int_t end_y_position;
+  Int_t value;
+};
+
+void AliL3HoughMaxFinder::FindAdaptedRowPeaks(Int_t kappawindow,Int_t xsize,Int_t ysize)
+{
+  
+  AliL3Histogram *hist = fCurrentHisto;
+  
+  if(!hist)
+    {
+      cerr<<"AliL3HoughMaxFinder : No histogram!"<<endl;
+      return;
+    }
+  
+  if(hist->GetNEntries() == 0)
+    return;
+  
+  Int_t xmin = hist->GetFirstXbin();
+  Int_t xmax = hist->GetLastXbin();
+  Int_t ymin = hist->GetFirstYbin();
+  Int_t ymax = hist->GetLastYbin();
+  
+  //Start by looking for pre-peaks:
+  
+  PreYPeak **local_maxima = new PreYPeak*[hist->GetNbinsY()];
+  
+  Short_t *nmaxs = new Short_t[hist->GetNbinsY()];
+  Int_t n2dmax = 0;
+  Int_t last_value,value;
+  for(Int_t ybin=ymin; ybin<=ymax; ybin++)
+    {
+      local_maxima[ybin-ymin] = new PreYPeak[hist->GetNbinsX()];
+      nmaxs[ybin-ymin] = 0;
+      last_value = 0;
+      Bool_t found = 0;
+      for(Int_t xbin=xmin; xbin<=xmax; xbin++)
+       {
+         value = hist->GetBinContent(hist->GetBin(xbin,ybin));
+         if(value >= fThreshold)
+           {
+             if(value > last_value)
+             //              if(value > (last_value + 1))
+               {
+                 local_maxima[ybin-ymin][nmaxs[ybin-ymin]].start_position = xbin;
+                 local_maxima[ybin-ymin][nmaxs[ybin-ymin]].end_position = xbin;
+                 local_maxima[ybin-ymin][nmaxs[ybin-ymin]].value = value;
+                 local_maxima[ybin-ymin][nmaxs[ybin-ymin]].prev_value = 0;
+                 found = 1;
+               }
+             if(value == last_value)
+               //            if(abs(value - last_value) <= 1)
+               if(found)
+                 local_maxima[ybin-ymin][nmaxs[ybin-ymin]].end_position = xbin;
+           }
+
+         if((value < fThreshold) || (value < last_value))
+           //    if((value < fThreshold) || (value < (last_value - 1)))
+           {
+             if(found)
+               {
+                 nmaxs[ybin-ymin]++;
+                 found = 0;
+               }
+           }
+         last_value = value;
+             
+       }
+      n2dmax += nmaxs[ybin-ymin];
+    }
+  
+  //  Pre2DPeak *maxima = new Pre2DPeak[n2dmax];
+  for(Int_t ybin=ymax; ybin >= ymin; ybin--)
+    {
+      for(Int_t i=0; i<nmaxs[ybin-ymin]; i++)
+       {
+         Int_t local_value = local_maxima[ybin-ymin][i].value;
+         Int_t local_prev_value = local_maxima[ybin-ymin][i].prev_value;
+         Int_t local_next_value = 0;
+
+         if(local_value<0)
+           continue; //already used
+
+         //start expanding in the psi-direction:
+
+         Int_t local_x_start = local_maxima[ybin-ymin][i].start_position;
+         Int_t local_x_end = local_maxima[ybin-ymin][i].end_position;
+         Int_t temp_x_start = local_maxima[ybin-ymin][i].start_position;
+         Int_t temp_x_end = local_maxima[ybin-ymin][i].end_position;
+
+         Int_t local_y=ybin-1,nybins=1;
+         
+         while(local_y >= ymin)
+           {
+             Bool_t found=0;
+             for(Int_t j=0; j<nmaxs[local_y-ymin]; j++)
+               {
+                 Int_t adapted_kappawindow;
+                 Float_t adapted_x,adapted_y;
+                 adapted_x = ((Float_t)local_x_start+(Float_t)local_x_end)/2.0;
+                 adapted_y = ((Float_t)ybin+(Float_t)(local_y+1))/2.0;
+                 //              if(adapted_y<(-0.46*adapted_x+86) && adapted_y>(-0.46*adapted_x+56) && adapted_y<(-0.215*adapted_x+68) && adapted_y>(-0.215*adapted_x+40))
+                 //              if(adapted_y<(-0.215*adapted_x+68) && adapted_y>(-0.215*adapted_x+40))
+                 //              if((adapted_y<(-0.46*adapted_x+86) && adapted_y>(-0.46*adapted_x+56)) && !(adapted_y<(-0.215*adapted_x+68) && adapted_y>(-0.215*adapted_x+40)))
+                 if(((adapted_y<(-0.46*adapted_x+86) && adapted_y>(-0.46*adapted_x+56)) && !(adapted_y<(-0.215*adapted_x+68) && adapted_y>(-0.215*adapted_x+40))) || 
+                    ((adapted_y<(-0.215*adapted_x+68) && adapted_y>(-0.215*adapted_x+40))&& !(adapted_y<(-0.46*adapted_x+86) && adapted_y>(-0.46*adapted_x+56))))
+                   adapted_kappawindow = 1;
+                 else
+                   adapted_kappawindow = 0;
+                 //temprorary here
+                 //              adapted_kappawindow = 0;
+                 
+                 if( (local_maxima[local_y-ymin][j].start_position <= (temp_x_end + kappawindow + adapted_kappawindow)) && (local_maxima[local_y-ymin][j].end_position >= temp_x_start)) 
+                   {
+                     if( local_maxima[local_y-ymin][j].value == local_value )
+                       {
+                         local_x_end = local_maxima[local_y-ymin][j].end_position;
+                         temp_x_start = local_maxima[local_y-ymin][j].start_position;
+                         temp_x_end = local_maxima[local_y-ymin][j].end_position;
+                         local_maxima[local_y-ymin][j].value = -1;
+                         found = 1;
+                         nybins++;
+                         break;
+                       }
+                     else
+                       {
+                         local_maxima[local_y-ymin][j].prev_value = local_value;
+                         local_next_value = local_maxima[local_y-ymin][j].value;
+                       }
+                   }
+               }
+             if(!found || local_y == ymin)//no more local maximas to be matched, so write the final peak and break the expansion:
+               {
+                 Int_t adapted_xsize;
+                 Float_t adapted_x,adapted_y;
+                 adapted_x = ((Float_t)local_x_start+(Float_t)local_x_end)/2.0;
+                 adapted_y = ((Float_t)ybin+(Float_t)(local_y+1))/2.0;
+                 //              if(adapted_y<(-0.46*adapted_x+86) && adapted_y>(-0.46*adapted_x+56) && adapted_y<(-0.215*adapted_x+68) && adapted_y>(-0.215*adapted_x+40))
+                 //              if(adapted_y<(-0.215*adapted_x+68) && adapted_y>(-0.215*adapted_x+40))
+                 if((adapted_y<(-0.46*adapted_x+86) && adapted_y>(-0.46*adapted_x+56)) && !(adapted_y<(-0.215*adapted_x+68) && adapted_y>(-0.215*adapted_x+40)))
+                   adapted_xsize = 1;
+                 else
+                   adapted_xsize = 2;
+                 //temprorary here
+                 //              adapted_xsize = 1;
+                 if((nybins > ysize) && ((local_x_end-local_x_start+1) > adapted_xsize) && (local_value > local_prev_value) && (local_value > local_next_value))
+                   {
+                     fXPeaks[fNPeaks] = hist->GetPreciseBinCenterX(adapted_x);
+                     fYPeaks[fNPeaks] = hist->GetPreciseBinCenterY(adapted_y);
+                     fWeight[fNPeaks] = local_value;
+#ifdef do_mc
+                     cout<<"Peak found at: "<<((Float_t)local_x_start+(Float_t)local_x_end)/2.0<<" "<<((Float_t)ybin+(Float_t)(local_y+1))/2.0<<" "<<fXPeaks[fNPeaks]<<" "<<fYPeaks[fNPeaks]<<" with weight "<<fWeight[fNPeaks]<<" and size "<<local_x_end-local_x_start+1<<" by "<<nybins<<endl;
+#endif
+                     fNPeaks++;
+                   }
+                 break;
+               }
+             else
+               local_y--;//Search continues...
+           }
+       }
+    }
+
+  for(Int_t i=0; i<hist->GetNbinsY(); i++)
+    delete local_maxima[i];
+
+  delete [] local_maxima;
+  delete [] nmaxs;
+  //delete [] maxima;
+}
 
 void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides)
 {
@@ -781,7 +967,6 @@ void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides)
     delete windowPt[i];
   delete [] windowPt;
   delete [] anotherPt;
-  
 }
 
 void AliL3HoughMaxFinder::SortPeaks(struct AxisWindow **a,Int_t first,Int_t last)
@@ -831,7 +1016,6 @@ Int_t AliL3HoughMaxFinder::PeakCompare(struct AxisWindow *a,struct AxisWindow *b
   if(a->weight < b->weight) return 1;
   if(a->weight > b->weight) return -1;
   return 0;
-
 }
 
 void AliL3HoughMaxFinder::FindPeak(Int_t t1,Double_t t2,Int_t t3)
@@ -1012,7 +1196,5 @@ void AliL3HoughMaxFinder::FindPeak(Int_t t1,Double_t t2,Int_t t3)
   delete [] m_up;
   
   //return track;
-
-    
 }
 
index fa78e96..39ce636 100644 (file)
@@ -53,6 +53,8 @@ class AliL3HoughMaxFinder {
   void FindBigMaxima();
   void FindMaxima(Int_t threshold=0);
   void FindAdaptedPeaks(Int_t nkappawindow,Float_t cut_ratio);
+  //Peak finder for HoughTransformerRow
+  void FindAdaptedRowPeaks(Int_t kappawindow,Int_t xsize,Int_t ysize);
   //More sophisticated peak finders:
   void FindPeak(Int_t t1,Double_t t2,Int_t t3);
   void FindPeak1(Int_t y_window=2,Int_t x_bin_sides=1);
index f61891b..2e9b785 100644 (file)
@@ -12,7 +12,7 @@
 #include "AliL3HoughMerger.h"
 #include "AliL3HoughTransformer.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 538e837..efda61a 100644 (file)
@@ -17,7 +17,7 @@
 #include <TH3.h>
 #include <TAxis.h>
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 3fac102..cea92d5 100644 (file)
@@ -9,14 +9,19 @@
 #include "AliL3HoughTrack.h"
 #include "AliL3Transform.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
+/** \class AliL3HoughTrack
+<pre>
 //_____________________________________________________________
 // AliL3HoughTrack
 //
 // Track class for Hough tracklets
+//
+</pre>
+*/
 
 ClassImp(AliL3HoughTrack)
 
@@ -34,14 +39,12 @@ AliL3HoughTrack::AliL3HoughTrack()
   fEta = 0;
 }
 
-
 AliL3HoughTrack::~AliL3HoughTrack()
 {
 }
 
 void AliL3HoughTrack::Set(AliL3Track *track)
 {
-  
   AliL3HoughTrack *tpt = (AliL3HoughTrack*)track;
   SetTrackParameters(tpt->GetKappa(),tpt->GetPsi(),tpt->GetWeight());
   SetEtaIndex(tpt->GetEtaIndex());
@@ -56,14 +59,17 @@ void AliL3HoughTrack::Set(AliL3Track *track)
   SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
   SetSlice(tpt->GetSlice());
   SetHits(tpt->GetNHits(),(UInt_t *)tpt->GetHitNumbers());
-  return;
-
-  fWeight = tpt->GetWeight();
-  fDLine = tpt->GetDLine();
-  fPsiLine = tpt->GetPsiLine();
-  SetNHits(tpt->GetWeight());
-  SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
-  fIsHelix = false;
+#ifdef ROWHOUGH
+  SetMCid(tpt->GetMCid());
+#endif
+  /*
+    fWeight = tpt->GetWeight();
+    fDLine = tpt->GetDLine();
+    fPsiLine = tpt->GetPsiLine();
+    SetNHits(tpt->GetWeight());
+    SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
+    fIsHelix = false;
+  */
 }
 
 Int_t AliL3HoughTrack::Compare(const AliL3Track *tpt) const
@@ -85,7 +91,6 @@ void AliL3HoughTrack::SetEta(Double_t f)
   SetTgl(tgl);
 }
 
-
 void AliL3HoughTrack::UpdateToFirstRow()
 {
   //Update the track parameters to the point where track cross
@@ -164,7 +169,6 @@ void AliL3HoughTrack::UpdateToFirstRow()
 
 void AliL3HoughTrack::SetTrackParameters(Double_t kappa,Double_t eangle,Int_t weight)
 {
-
   fWeight = weight;
   fMinDist = 100000;
   SetKappa(kappa);
@@ -202,23 +206,19 @@ void AliL3HoughTrack::SetLineParameters(Double_t psi,Double_t D,Int_t weight,Int
   SetNHits(1);
   SetRowRange(rowrange[0],rowrange[1]);
   fIsHelix = false;
-
 }
 
 void AliL3HoughTrack::SetBestMCid(Int_t mcid,Double_t min_dist)
 {
-  
   if(min_dist < fMinDist)
     {
       fMinDist = min_dist;
       SetMCid(mcid);
     }
-  
 }
 
 void AliL3HoughTrack::GetLineCrossingPoint(Int_t padrow,Float_t *xy)
 {
-  
   if(fIsHelix)
     {
       printf("AliL3HoughTrack::GetLineCrossingPoint : Track is not a line\n");
@@ -231,5 +231,4 @@ void AliL3HoughTrack::GetLineCrossingPoint(Int_t padrow,Float_t *xy)
   Float_t yhit = a*xhit + b;
   xy[0] = xhit;
   xy[1] = yhit;
-
 }
index b501875..19e3afb 100644 (file)
@@ -12,7 +12,7 @@
 #include "AliL3DigitData.h"
 #include "AliL3HistogramAdaptive.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
diff --git a/HLT/hough/AliL3HoughTransformerGap.cxx b/HLT/hough/AliL3HoughTransformerGap.cxx
deleted file mode 100644 (file)
index 0f78ace..0000000
+++ /dev/null
@@ -1,804 +0,0 @@
-// @(#) $Id$
-
-// Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
-//*-- Copyright &copy ALICE HLT Group
-
-#include "AliL3StandardIncludes.h"
-
-#include "AliL3Logging.h"
-#include "AliL3MemHandler.h"
-#include "AliL3Transform.h"
-#include "AliL3DigitData.h"
-#include "AliL3HistogramAdaptive.h"
-#include "AliL3HoughTransformerGap.h"
-
-#if GCCVERSION == 3
-using namespace std;
-#endif
-
-//_____________________________________________________________
-// AliL3HoughTransformerGap
-//
-// Hough transformation class
-//
-
-ClassImp(AliL3HoughTransformerGap)
-
-#ifdef COUNTINGROWS
-UChar_t **AliL3HoughTransformerGap::fRowCount = 0;
-UChar_t **AliL3HoughTransformerGap::fGapCount = 0;
-UChar_t **AliL3HoughTransformerGap::fCurrentRowCount = 0;
-#endif
-
-AliL3HoughTransformerGap::AliL3HoughTransformerGap()
-{
-  //Default constructor
-  fParamSpace = 0;
-  fDoMC = kFALSE;;
-
-#ifdef do_mc
-  fTrackID = 0;
-#endif
-  fZVertex = 0.;
-#ifdef SUBSLICES
-  fEtaSpace = 0;
-  fEtaSpaceSize = 0;
-#endif
-#ifdef COUNTINGROWS
-  fLUT2sinphi0up=0;
-  fLUT2sinphi0low=0;
-  fLUT2cosphi0up=0;
-  fLUT2cosphi0low=0;
-#endif
-}
-
-AliL3HoughTransformerGap::AliL3HoughTransformerGap(Int_t slice,Int_t patch,Int_t n_eta_segments,Bool_t DoEtaOverlap,Bool_t DoMC,Float_t zvertex) : AliL3HoughBaseTransformer(slice,patch,n_eta_segments)
-{
-  //Normal constructor
-  fParamSpace = 0;
-  fDoMC = kFALSE;
-  fDoMC=kFALSE;
-
-#ifdef do_mc
-  fTrackID = 0;
-#endif
-  fZVertex = zvertex;
-#ifdef SUBSLICES
-  fEtaSpace = 0;
-  fEtaSpaceSize = 800/n_eta_segments;
-#endif
-#ifdef COUNTINGROWS
-  fLUT2sinphi0up=0;
-  fLUT2sinphi0low=0;
-  fLUT2cosphi0up=0;
-  fLUT2cosphi0low=0;
-#endif
-}
-
-AliL3HoughTransformerGap::~AliL3HoughTransformerGap()
-{
-  DeleteHistograms();
-#ifdef do_mc
-  if(fTrackID)
-    {
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       {
-         if(!fTrackID[i]) continue;
-         delete fTrackID[i];
-       }
-      delete [] fTrackID;
-    }
-#endif
-#ifdef SUBSLICES
-  if(fEtaSpace)
-    {
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       {
-         if(!fEtaSpace[i]) continue;
-         delete fEtaSpace[i];
-       }
-      delete [] fEtaSpace;
-    }
-#endif
-#ifdef COUNTINGROWS
-  if(fRowCount)
-    {
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       {
-         if(!fRowCount[i]) continue;
-         delete fRowCount[i];
-       }
-      delete [] fRowCount;
-    }
-  if(fGapCount)
-    {
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       {
-         if(!fGapCount[i]) continue;
-         delete fGapCount[i];
-       }
-      delete [] fGapCount;
-    }
-  if(fCurrentRowCount)
-    {
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       {
-         if(fCurrentRowCount[i])
-           delete fCurrentRowCount[i];
-       }
-      delete [] fCurrentRowCount;
-    }
-#endif
-}
-
-void AliL3HoughTransformerGap::DeleteHistograms()
-{
-#ifdef COUNTINGROWS
-  delete[] fLUT2sinphi0up;
-  delete[] fLUT2cosphi0up;
-  delete[] fLUT2sinphi0low;
-  delete[] fLUT2cosphi0low;
-
-  fLUT2sinphi0up=0;
-  fLUT2cosphi0up=0;
-  fLUT2sinphi0low=0;
-  fLUT2cosphi0low=0;
-#endif
-  if(!fParamSpace)
-    return;
-  for(Int_t i=0; i<GetNEtaSegments(); i++)
-    {
-      if(!fParamSpace[i]) continue;
-      delete fParamSpace[i];
-    }
-  delete [] fParamSpace;
-}
-
-void AliL3HoughTransformerGap::CreateHistograms(Int_t nxbin,Float_t pt_min,
-                                            Int_t nybin,Float_t phimin,Float_t phimax)
-{
-  //Create the histograms (parameter space).
-  //These are 2D histograms, span by kappa (curvature of track) 
-  //and phi0 (emission angle with x-axis).
-  //The arguments give the range and binning; 
-  //nxbin = #bins in kappa
-  //nybin = #bins in phi0
-  //pt_min = mimium Pt of track (corresponding to maximum kappa)
-  //phi_min = mimimum phi0 
-  //phi_max = maximum phi0 
-    
-  Double_t x = AliL3Transform::GetBFact()*AliL3Transform::GetBField()/pt_min;
-  //Double_t torad = AliL3Transform::Pi()/180;
-  CreateHistograms(nxbin,-1.*x,x,nybin,phimin/**torad*/,phimax/**torad*/);
-}
-
-void AliL3HoughTransformerGap::CreateHistograms(Int_t nxbin,Float_t xmin,Float_t xmax,
-                                            Int_t nybin,Float_t ymin,Float_t ymax)
-{
-  
-  fParamSpace = new AliL3Histogram*[GetNEtaSegments()];
-  
-  Char_t histname[256];
-  for(Int_t i=0; i<GetNEtaSegments(); i++)
-    {
-      sprintf(histname,"paramspace_%d",i);
-      fParamSpace[i] = new AliL3Histogram(histname,"",nxbin,xmin,xmax,nybin,ymin,ymax);
-    }
-  
-#ifdef do_mc
-  if(fDoMC)
-    {
-      AliL3Histogram *hist = fParamSpace[0];
-      Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
-      cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(TrackIndex)<<" bytes to fTrackID"<<endl;
-      fTrackID = new TrackIndex*[GetNEtaSegments()];
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       fTrackID[i] = new TrackIndex[ncells];
-    }
-#endif
-#ifdef SUBSLICES
-  AliL3Histogram *hist = fParamSpace[0];
-  Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
-  cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*fEtaSpaceSize*sizeof(Int_t)<<" bytes to fEtaSpace"<<endl;
-  fEtaSpace = new Int_t*[GetNEtaSegments()];
-  for(Int_t i=0; i<GetNEtaSegments(); i++)
-    fEtaSpace[i] = new Int_t[ncells*fEtaSpaceSize];
-#endif
-#ifdef COUNTINGROWS
-  AliL3Histogram *hist = fParamSpace[0];
-  Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
-  if(!fRowCount)
-    {
-      cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fRowCount"<<endl;
-      fRowCount = new UChar_t*[GetNEtaSegments()];
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       fRowCount[i] = new UChar_t[ncells];
-    }
-  if(!fGapCount)
-    {
-      cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fGapCount"<<endl;
-      fGapCount = new UChar_t*[GetNEtaSegments()];
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       fGapCount[i] = new UChar_t[ncells];
-    }
-  if(!fCurrentRowCount)
-    {
-      cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fCurrentRowCount"<<endl;
-      fCurrentRowCount = new UChar_t*[GetNEtaSegments()];
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       fCurrentRowCount[i] = new UChar_t[ncells];
-    }
-
-  //create lookup table for sin and cos
-  Int_t nbins = hist->GetNbinsY()+2;
-  fLUT2sinphi0up=new Float_t[nbins];
-  fLUT2cosphi0up=new Float_t[nbins];
-  fLUT2sinphi0low=new Float_t[nbins];
-  fLUT2cosphi0low=new Float_t[nbins];
-  Float_t hist_bin=hist->GetBinWidthY()/2.0;
-  for(Int_t i=hist->GetFirstYbin(); i<=hist->GetLastYbin(); i++){
-    Float_t phi0=hist->GetBinCenterY(i);
-    fLUT2sinphi0low[i]=2.*sin(phi0+hist_bin);
-    fLUT2cosphi0low[i]=2.*cos(phi0+hist_bin);
-    fLUT2sinphi0up[i]=2.*sin(phi0-hist_bin);
-    fLUT2cosphi0up[i]=2.*cos(phi0-hist_bin);
-  }  
-#endif
-
-}
-
-void AliL3HoughTransformerGap::Reset()
-{
-  //Reset all the histograms
-
-  if(!fParamSpace)
-    {
-      LOG(AliL3Log::kWarning,"AliL3HoughTransformer::Reset","Histograms")
-       <<"No histograms to reset"<<ENDLOG;
-      return;
-    }
-  
-  for(Int_t i=0; i<GetNEtaSegments(); i++)
-    fParamSpace[i]->Reset();
-  
-#ifdef do_mc
-  if(fDoMC)
-    {
-      AliL3Histogram *hist = fParamSpace[0];
-      Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
-      for(Int_t i=0; i<GetNEtaSegments(); i++)
-       memset(fTrackID[i],0,ncells*sizeof(TrackIndex));
-    }
-#endif
-#ifdef SUBSLICES
-  AliL3Histogram *hist = fParamSpace[0];
-  Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
-  for(Int_t i=0; i<GetNEtaSegments(); i++)
-    memset(fEtaSpace[i],0,ncells*fEtaSpaceSize*sizeof(Int_t));
-#endif
-#ifdef COUNTINGROWS
-  AliL3Histogram *hist = fParamSpace[0];
-  Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
-  for(Int_t i=0; i<GetNEtaSegments(); i++)
-    {
-      memset(fRowCount[i],0,ncells*sizeof(UChar_t));
-      memset(fGapCount[i],1,ncells*sizeof(UChar_t));
-      memset(fCurrentRowCount[i],160,ncells*sizeof(UChar_t));
-    }
-#endif
-}
-
-Int_t AliL3HoughTransformerGap::GetEtaIndex(Double_t eta)
-{
-  //Return the histogram index of the corresponding eta. 
-
-  Double_t etaslice = (GetEtaMax() - GetEtaMin())/GetNEtaSegments();
-  Double_t index = (eta-GetEtaMin())/etaslice;
-  return (Int_t)index;
-}
-
-inline AliL3Histogram *AliL3HoughTransformerGap::GetHistogram(Int_t eta_index)
-{
-  if(!fParamSpace || eta_index >= GetNEtaSegments() || eta_index < 0)
-    return 0;
-  if(!fParamSpace[eta_index])
-    return 0;
-  return fParamSpace[eta_index];
-}
-
-Double_t AliL3HoughTransformerGap::GetEta(Int_t eta_index,Int_t slice)
-{
-  Double_t eta_slice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments();
-  Double_t eta=0;
-  eta=(Double_t)((eta_index+0.5)*eta_slice);
-  return eta;
-}
-
-#ifdef SUBSLICES
-Int_t AliL3HoughTransformerGap::GetEtaSubIndex(Double_t eta,Int_t eta_index)
-{
-  //Return the eta space sub index of the corresponding eta. 
-
-  Double_t sub_eta = GetEta(eta_index,GetSlice());
-  Int_t eta_subindex = (Int_t)((eta-sub_eta+0.5/(Double_t)GetNEtaSegments())*800.);
-  return eta_subindex;
-}
-
-Double_t AliL3HoughTransformerGap::GetMaxSubEta(Float_t kappa,Float_t phi0,Int_t eta_index,Int_t slice)
-{
-  AliL3Histogram *hist = fParamSpace[eta_index];
-  Int_t bin = hist->FindBin(kappa,phi0);
-  Int_t maxbin=fEtaSpaceSize/2;
-  Int_t maxvalue=0;
-  for(Int_t i=0;i<fEtaSpaceSize;i++) {
-    if(fEtaSpace[eta_index][bin*fEtaSpaceSize+i]>maxvalue) {
-      maxbin=i;
-      maxvalue=fEtaSpace[eta_index][bin*fEtaSpaceSize+i];
-    }
-  }
-  Double_t eta = GetEta(eta_index,slice)+((Double_t)maxbin-0.5*(Double_t)fEtaSpaceSize+0.5)/800.;
-  cout<<"DEBUG maxetasubindex "<<eta<<" "<<eta_index<<" "<<maxbin<<endl;
-  return eta;
-}
-#endif
-
-#ifdef notused
-void AliL3HoughTransformerGap::TransformCircle()
-{
-
-  AliL3DigitRowData *tempPt = GetDataPointer();
-  if(!tempPt)
-    {
-      LOG(AliL3Log::kError,"AliL3HoughTransformer::TransformCircle","Data")
-       <<"No input data "<<ENDLOG;
-      return;
-    }
-
-  Int_t ipatch = GetPatch();
-  //Loop over the padrows:
-  for(Int_t i=AliL3Transform::GetFirstRow(GetPatch()); i<=AliL3Transform::GetLastRow(GetPatch()); i++)
-    {
-      //Get the data on this padrow:
-      AliL3DigitData *digPt = tempPt->fDigitData;
-      if(i != (Int_t)tempPt->fRow)
-       {
-         cerr<<"AliL3HoughTransform::TransformCircle : Mismatching padrow numbering "<<i<<" "<<(Int_t)tempPt->fRow<<endl;
-         continue;
-       }
-      
-      //Loop over the data on this padrow:
-      for(UInt_t j=0; j<tempPt->fNDigit; j++)
-       {
-         UShort_t charge = digPt[j].fCharge;
-         UChar_t pad = digPt[j].fPad;
-         UShort_t time = digPt[j].fTime;
-         if((Int_t)charge <= GetLowerThreshold())
-           continue;
-         
-         if((Int_t)charge > GetUpperThreshold())
-           charge = GetUpperThreshold();
-         
-         Int_t sector,row;
-         Float_t xyz[3];
-
-         //Transform data to local cartesian coordinates:
-         AliL3Transform::Slice2Sector(GetSlice(),i,sector,row);
-         AliL3Transform::Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
-                 
-         //Calculate the eta:
-         xyz[2] -= fZVertex;
-         Double_t eta = AliL3Transform::GetEta(xyz);
-         
-         //Get the corresponding index, which determines which histogram to fill:
-         Int_t eta_index = GetEtaIndex(eta);
-                 
-         if(eta_index < 0 || eta_index >= GetNEtaSegments())
-           continue;
-         
-         //Get the correct histogrampointer:
-         AliL3Histogram *hist = fParamSpace[eta_index];
-         if(!hist)
-           {
-             cerr<<"AliL3HoughTransformer::TransformCircle : Error getting histogram in index "<<eta_index<<endl;
-             continue;
-           }
-         
-         //Do the transformation:
-#ifndef COUNTINGROWS
-         Float_t R = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]); 
-         Float_t phi = AliL3Transform::GetPhi(xyz);
-#else
-         xyz[1] -= AliL3Transform::GetPadPitchWidth(ipatch)/2.0;
-         Float_t R1 = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]); 
-         Float_t phi1 = AliL3Transform::GetPhi(xyz);
-         xyz[1] += AliL3Transform::GetPadPitchWidth(ipatch);
-         Float_t R2 = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]); 
-         Float_t phi2 = AliL3Transform::GetPhi(xyz);
-#endif
-#ifdef SUBSLICES         
-         Int_t eta_subindex = GetEtaSubIndex(eta,eta_index);
-#endif   
-         //Fill the histogram along the phirange
-         for(Int_t b=hist->GetFirstYbin(); b<=hist->GetLastYbin(); b++)
-           {
-             Float_t phi0 = hist->GetBinCenterY(b);
-#ifndef COUNTINGROWS
-             Float_t kappa = 2*sin(phi - phi0)/R;
-#else
-             Float_t kappa1 = 2*sin(phi1 - phi0 - hist->GetBinWidthY()/2.)/R1;
-             Float_t kappa2 = 2*sin(phi2 - phi0 + hist->GetBinWidthY()/2.)/R2;
-#endif
-             //hist->Fill(kappa,phi0,(int)rint(log((Float_t)charge)));
-             //              hist->Fill(kappa,phi0,charge);
-#ifndef COUNTINGROWS
-             hist->Fill(kappa,phi0,30);
-             //hist->Fill(kappa,phi0,1);
-#else
-             Int_t binx1,binx2;
-             binx1 = hist->FindXbin(kappa1);
-             binx2 = hist->FindXbin(kappa2);
-
-             //              cout<<"DEBUG counting rows "<<ipatch<<" "<<i<<" "<<" bins "<<binx1<<" "<<binx2<<" "<<" kappas "<<kappa1<<" "<<kappa2<<endl;
-
-             Int_t bin1 = hist->GetBin(binx1,b);
-             Int_t bin2 = hist->GetBin(binx2,b);
-             if((bin1==0) || (bin2==0)) continue;
-             //              cout<<"DEBUG counting rows "<<ipatch<<" "<<i<<" "<<" bins "<<binx1<<" "<<binx2<<" "<<bin1<<" "<<bin2<<" kappas "<<kappa1<<" "<<kappa2<<endl;
-             for(Int_t bin=bin1;bin<=bin2;bin++) {
-               fRowCount[eta_index][bin]=fRowCount[eta_index][bin]|(1<<(i-AliL3Transform::GetFirstRow(ipatch)));
-             }
-#endif
-#ifndef COUNTINGROWS
-#ifdef SUBSLICES
-             //Fill Eta space histogram
-             Int_t bin = hist->FindBin(kappa,phi0);
-             //              cout<<" DEBUG etasubindex "<<eta<<" "<<eta_index<<" "<<eta_subindex<<endl;
-             fEtaSpace[eta_index][bin*fEtaSpaceSize+eta_subindex]++;
-#endif
-#ifdef do_mc
-             if(fDoMC)
-               {
-                 Int_t bin = hist->FindBin(kappa,phi0);
-                 for(Int_t t=0; t<3; t++)
-                   {
-                     Int_t label = digPt[j].fTrackID[t];
-                     if(label < 0) break;
-                     UInt_t c;
-                     for(c=0; c<MaxTrack; c++)
-                       if(fTrackID[eta_index][bin].fLabel[c] == label || fTrackID[eta_index][bin].fNHits[c] == 0)
-                         break;
-                     if(c == MaxTrack-1) cerr<<"AliL3HoughTransformer::TransformCircle : Array reached maximum!! "<<c<<endl;
-                     fTrackID[eta_index][bin].fLabel[c] = label;
-                     fTrackID[eta_index][bin].fNHits[c]++;
-                   }
-               }
-#endif
-#endif
-           }
-       }
-      
-      //Move the data pointer to the next padrow:
-      AliL3MemHandler::UpdateRowPointer(tempPt);
-    }
-}
-#endif
-
-#ifdef COUNTINGROWS
-
-struct EtaRow {
-  UChar_t start_pad;
-  UChar_t end_pad;
-  Bool_t found;
-  Float_t start_y;
-};
-
-void AliL3HoughTransformerGap::TransformCircle()
-{
-  //Assumes that all the etaslice histos are the same!
-  AliL3Histogram *h = fParamSpace[0];
-  //Float_t hist_bin = h->GetBinWidthY()/2.;
-  Int_t first_bin = h->GetFirstYbin();
-  Int_t last_bin = h->GetLastYbin();
-  Float_t x_min = h->GetXmin();
-  Float_t x_max = h->GetXmax();
-  Float_t x_bin = (x_max-x_min)/h->GetNbinsX();
-  Int_t first_binx = h->GetFirstXbin()-1;
-  Int_t last_binx = h->GetLastXbin()+1;
-  Int_t nbinx = h->GetNbinsX()+2;
-
-  UChar_t last_pad;
-  Int_t last_eta_index=-1;;
-  EtaRow *eta_clust = new EtaRow[GetNEtaSegments()];
-
-  AliL3DigitRowData *tempPt = GetDataPointer();
-  if(!tempPt)
-    {
-      LOG(AliL3Log::kError,"AliL3HoughTransformer::TransformCircle","Data")
-       <<"No input data "<<ENDLOG;
-      return;
-    }
-
-  Int_t ipatch = GetPatch();
-  Double_t pad_pitch = AliL3Transform::GetPadPitchWidth(ipatch);
-  Int_t islice = GetSlice();
-  //Loop over the padrows:
-  for(Int_t i=AliL3Transform::GetFirstRow(ipatch); i<=AliL3Transform::GetLastRow(ipatch); i++)
-    {
-      Int_t npads = AliL3Transform::GetNPads(i)-1;
-
-      last_pad = 255;
-      //Flush eta clusters array
-      memset(eta_clust,0,GetNEtaSegments()*sizeof(EtaRow));  
-
-      Float_t x = AliL3Transform::Row2X(i);
-      Float_t x2 = x*x;
-      Float_t y=0.,r2=0.;
-
-      UChar_t lrow = i-AliL3Transform::GetFirstRow(ipatch);
-
-      //Get the data on this padrow:
-      AliL3DigitData *digPt = tempPt->fDigitData;
-      if(i != (Int_t)tempPt->fRow)
-       {
-         cerr<<"AliL3HoughTransform::TransformCircle : Mismatching padrow numbering "<<i<<" "<<(Int_t)tempPt->fRow<<endl;
-         continue;
-       }
-      //      cout<<" Starting row "<<i<<endl;
-      //Loop over the data on this padrow:
-      for(UInt_t j=0; j<tempPt->fNDigit; j++)
-       {
-         UShort_t charge = digPt[j].fCharge;
-         if((Int_t)charge <= GetLowerThreshold())
-           continue;
-         UChar_t pad = digPt[j].fPad;
-         UShort_t time = digPt[j].fTime;
-         if(pad != last_pad)
-           {
-             y = (pad-0.5*npads)*pad_pitch;
-             Float_t y2 = y*y;
-             r2 = x2 + y2;
-             last_eta_index = -1;
-           }
-
-         //Transform data to local cartesian coordinates:
-         Float_t z = AliL3Transform::GetZFast(islice,(Int_t)time,fZVertex);
-         Float_t z2 = z*z;
-         //Calculate the eta:
-         Double_t r = sqrt(r2+z2);
-         Double_t eta = 0.5 * log((r+z)/(r-z));
-         //Get the corresponding index, which determines which histogram to fill:
-         Int_t eta_index = GetEtaIndex(eta);
-
-         if(eta_index == last_eta_index) continue;
-
-         //      cout<<" Digit at patch "<<ipatch<<" row "<<i<<" pad "<<(Int_t)pad<<" time "<<time<<" etaslice "<<eta_index<<endl;
-         
-         if(eta_index < 0 || eta_index >= GetNEtaSegments())
-           continue;
-
-         if(!eta_clust[eta_index].found)
-           {
-             eta_clust[eta_index].start_pad = pad;
-             eta_clust[eta_index].end_pad = pad;
-             eta_clust[eta_index].start_y = y - pad_pitch/2.0;
-             eta_clust[eta_index].found = 1;
-             continue;
-           }
-         else
-           {
-             if(pad <= (eta_clust[eta_index].end_pad + 1))
-               eta_clust[eta_index].end_pad = pad;
-             else
-               {
-
-                 //              cout<<" Cluster found at etaslice "<<eta_index<<" from pad "<<(Int_t)eta_clust[eta_index].start_pad<<" to pad "<<(Int_t)eta_clust[eta_index].end_pad<<endl;
-                 //Get the correct histogrampointer:
-                 AliL3Histogram *hist = fParamSpace[eta_index];
-                 if(!hist)
-                   {
-                     cerr<<"AliL3HoughTransformer::TransformCircle : Error getting histogram in index "<<eta_index<<endl;
-                     continue;
-                   }
-                 UChar_t *nrows = fRowCount[eta_index];
-                 UChar_t *ngaps = fGapCount[eta_index];
-                 UChar_t *currentrow = fCurrentRowCount[eta_index];
-
-                 //Do the transformation:
-                 Float_t start_y = eta_clust[eta_index].start_y;
-                 Float_t R1 = x2 + start_y*start_y; 
-                 Float_t end_y = (eta_clust[eta_index].end_pad-0.5*npads)*pad_pitch + pad_pitch/2.0;
-                 Float_t R2 = x2 + end_y*end_y; 
-
-                 //Fill the histogram along the phirange
-                 for(Int_t b=first_bin; b<=last_bin; b++)
-                   {
-#ifdef WOLUT
-                     Float_t phi0 = hist->GetBinCenterY(b);
-                     //                      Float_t kappa1 = 2*sin(phi1 - phi0 - hist_bin)/R1;
-                     //                      Float_t kappa2 = 2*sin(phi2 - phi0 + hist_bin)/R2;
-                     Float_t kappa1 = 2*(start_y*cos(phi0 + hist_bin)-x*sin(phi0 + hist_bin))/R1;
-                     if(kappa1>x_max) continue;
-                     Float_t kappa2 = 2*(end_y*cos(phi0 - hist_bin)-x*sin(phi0 - hist_bin))/R2;
-                     if(kappa2<x_min) break;
-#else
-                     Float_t kappa1 = (start_y*fLUT2cosphi0low[b]-x*fLUT2sinphi0low[b])/R1;
-                     if(kappa1>x_max) continue;
-                     Float_t kappa2 = (end_y*fLUT2cosphi0up[b]-x*fLUT2sinphi0up[b])/R2;
-                     if(kappa2<x_min) break;
-#endif
-                     Int_t binx1,binx2;
-                     if(kappa1<x_min)
-                       binx1 = first_binx;
-                     else
-                       binx1 = 1 + (Int_t)((kappa1-x_min)/x_bin);
-                     if(kappa2>x_max)
-                       binx2 = last_binx;
-                     else
-                       binx2 = 1 + (Int_t)((kappa2-x_min)/x_bin);
-                     
-                     Int_t temp_bin = b*nbinx + binx1;
-                     UChar_t *nrows2 = nrows + temp_bin;
-                     UChar_t *ngaps2 = ngaps + temp_bin;
-                     UChar_t *currentrow2 = currentrow + temp_bin;
-                     for(Int_t bin=binx1;bin<=binx2;bin++)
-                       {
-                         //Assumes threshold about 32
-                         if(*ngaps2 < 5) {
-                           if(lrow != *currentrow2)
-                             {
-                               (*nrows2)++;
-                               if(lrow > (*currentrow2 + 1))
-                                 (*ngaps2)++;
-                               *currentrow2=lrow;
-                             }
-                         }
-                         nrows2++;
-                         ngaps2++;
-                         currentrow2++;
-                       }
-                   }
-                 //End of the transformation
-
-                 eta_clust[eta_index].start_pad = pad;
-                 eta_clust[eta_index].end_pad = pad;
-                 eta_clust[eta_index].start_y = y - pad_pitch/2.0;
-               }
-           }
-         last_pad = pad;
-         last_eta_index = eta_index;
-       }
-      //Write remaining clusters
-      for(Int_t eta_index = 0;eta_index < GetNEtaSegments();eta_index++)
-       {
-         //Check for empty row
-         if((eta_clust[eta_index].start_pad == 0) && (eta_clust[eta_index].end_pad == 0)) continue;
-
-         //Get the correct histogrampointer:
-         AliL3Histogram *hist = fParamSpace[eta_index];
-         if(!hist)
-           {
-             cerr<<"AliL3HoughTransformer::TransformCircle : Error getting histogram in index "<<eta_index<<endl;
-             continue;
-           }
-         UChar_t *nrows = fRowCount[eta_index];
-         UChar_t *ngaps = fGapCount[eta_index];
-         UChar_t *currentrow = fCurrentRowCount[eta_index];
-
-         //Do the transformation:
-         Float_t start_y = eta_clust[eta_index].start_y;
-         Float_t R1 = x2 + start_y*start_y; 
-         Float_t end_y = (eta_clust[eta_index].end_pad-0.5*npads)*pad_pitch + pad_pitch/2.0;
-         Float_t R2 = x2 + end_y*end_y; 
-
-         //Fill the histogram along the phirange
-         for(Int_t b=first_bin; b<=last_bin; b++)
-           {
-#ifdef WOLUT
-             Float_t phi0 = hist->GetBinCenterY(b);
-             //                      Float_t kappa1 = 2*sin(phi1 - phi0 - hist_bin)/R1;
-             //                      Float_t kappa2 = 2*sin(phi2 - phi0 + hist_bin)/R2;
-             Float_t kappa1 = 2*(start_y*cos(phi0 + hist_bin)-x*sin(phi0 + hist_bin))/R1;
-             if(kappa1>x_max) continue;
-             Float_t kappa2 = 2*(end_y*cos(phi0 - hist_bin)-x*sin(phi0 - hist_bin))/R2;
-             if(kappa2<x_min) break;
-#else
-             Float_t kappa1 = (start_y*fLUT2cosphi0low[b]-x*fLUT2sinphi0low[b])/R1;
-             if(kappa1>x_max) continue;
-             Float_t kappa2 = (end_y*fLUT2cosphi0up[b]-x*fLUT2sinphi0up[b])/R2;
-             if(kappa2<x_min) break;
-#endif
-             Int_t binx1,binx2;
-             if(kappa1<x_min)
-               binx1 = first_binx;
-             else
-               binx1 = 1 + (Int_t)((kappa1-x_min)/x_bin);
-             if(kappa2>x_max)
-               binx2 = last_binx;
-             else
-               binx2 = 1 + (Int_t)((kappa2-x_min)/x_bin);
-
-             Int_t temp_bin = b*nbinx + binx1;
-             UChar_t *nrows2 = nrows + temp_bin;
-             UChar_t *ngaps2 = ngaps + temp_bin;
-             UChar_t *currentrow2 = currentrow + temp_bin;
-             for(Int_t bin=binx1;bin<=binx2;bin++)
-               {
-                 //Assumes threshold about 32
-                 if(*ngaps2 < 5) {
-                   if(lrow != *currentrow2)
-                     {
-                       (*nrows2)++;
-                       if(lrow > (*currentrow2 + 1))
-                         (*ngaps2)++;
-                       *currentrow2=lrow;
-                     }
-                 }
-                 nrows2++;
-                 ngaps2++;
-                 currentrow2++;
-               }
-           }
-         //End of the transformation
-       }
-
-      //Move the data pointer to the next padrow:
-      AliL3MemHandler::UpdateRowPointer(tempPt);
-    }
-  delete [] eta_clust;
-}
-#endif
-
-#ifdef COUNTINGROWS
-/*
-void AliL3HoughTransformer::CalculateNRows()
-{
-  Int_t ipatch = GetPatch();
-
-  for(Int_t i=0; i<GetNEtaSegments(); i++) {
-    AliL3Histogram *hist = fParamSpace[i];
-
-    Int_t xmin = hist->GetFirstXbin();
-    Int_t xmax = hist->GetLastXbin();
-    Int_t ymin = hist->GetFirstYbin();
-    Int_t ymax = hist->GetLastYbin();
-    Int_t sum,ngaps;
-    Bool_t isPrevRow;
-
-    for(Int_t ybin=ymin; ybin<=ymax; ybin++)
-      {
-      for(Int_t xbin=xmin; xbin<=xmax; xbin++)
-       {
-         sum = 0;
-         ngaps = 1;
-         isPrevRow = kFALSE;
-         Int_t bin = hist->GetBin(xbin,ybin);
-         for(Int_t irow=0;irow<=(AliL3Transform::GetLastRow(ipatch)-AliL3Transform::GetFirstRow(ipatch));irow++) {
-           if(((fRowCount[i][bin]>>irow)&1)==1) {
-             sum++;
-             isPrevRow = kTRUE;
-           }
-           else {
-             if(isPrevRow)
-               ngaps++;
-             isPrevRow = kFALSE;
-           }
-         }
-         hist->SetBinContent(bin,sum/ngaps);
-       }
-      }
-    hist->SetNEntries(1);
-  }
-}
-*/
-
-UChar_t *AliL3HoughTransformerGap::GetRowCount(Int_t eta_index)
-{
-  return fRowCount[eta_index];
-}
-
-UChar_t *AliL3HoughTransformerGap::GetGapCount(Int_t eta_index)
-{
-  return fGapCount[eta_index];
-}
-
-#endif
index 79ddabc..1eaadf2 100644 (file)
@@ -19,7 +19,7 @@
 
 #include <TH2.h>
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index d9129d7..d458201 100644 (file)
@@ -13,7 +13,7 @@
 #include "AliL3Histogram.h"
 #include "AliL3HistogramAdaptive.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index fa2904b..24c37cd 100644 (file)
@@ -17,7 +17,7 @@
 #include "AliL3Vertex.h"
 #include "AliL3Fitter.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
diff --git a/HLT/hough/AliL3HoughTransformerRow.cxx b/HLT/hough/AliL3HoughTransformerRow.cxx
new file mode 100644 (file)
index 0000000..80c48c0
--- /dev/null
@@ -0,0 +1,776 @@
+// @(#) $Id$
+
+// Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
+//*-- Copyright &copy ALICE HLT Group
+
+#include "AliL3StandardIncludes.h"
+
+#include "AliL3Logging.h"
+#include "AliL3MemHandler.h"
+#include "AliL3Transform.h"
+#include "AliL3DigitData.h"
+#include "AliL3HistogramAdaptive.h"
+#include "AliL3HoughTransformerRow.h"
+
+#if __GNUC__ == 3
+using namespace std;
+#endif
+
+//_____________________________________________________________
+// AliL3HoughTransformerRow
+//
+// TPC rows Hough transformation class
+//
+
+ClassImp(AliL3HoughTransformerRow)
+
+UChar_t **AliL3HoughTransformerRow::fRowCount = 0;
+UChar_t **AliL3HoughTransformerRow::fGapCount = 0;
+UChar_t **AliL3HoughTransformerRow::fCurrentRowCount = 0;
+#ifdef do_mc
+TrackIndex **AliL3HoughTransformerRow::fTrackID = 0;
+#endif
+
+
+AliL3HoughTransformerRow::AliL3HoughTransformerRow()
+{
+  //Default constructor
+  fParamSpace = 0;
+  fDoMC = kFALSE;;
+
+  fLUT2sinphi0up=0;
+  fLUT2sinphi0low=0;
+  fLUT2cosphi0up=0;
+  fLUT2cosphi0low=0;
+
+  fLUTforwardZ=0;
+  fLUTforwardZ2=0;
+  fLUTbackwardZ=0;
+  fLUTbackwardZ2=0;
+}
+
+AliL3HoughTransformerRow::AliL3HoughTransformerRow(Int_t slice,Int_t patch,Int_t n_eta_segments,Bool_t DoMC,Float_t zvertex) : AliL3HoughBaseTransformer(slice,patch,n_eta_segments,zvertex)
+{
+  //Normal constructor
+  fParamSpace = 0;
+  fDoMC = kFALSE;
+  fDoMC=kFALSE;
+#ifdef do_mc
+  fDoMC = kTRUE;
+#endif
+
+  fLUT2sinphi0up=0;
+  fLUT2sinphi0low=0;
+  fLUT2cosphi0up=0;
+  fLUT2cosphi0low=0;
+
+  fLUTforwardZ=0;
+  fLUTforwardZ2=0;
+  fLUTbackwardZ=0;
+  fLUTbackwardZ2=0;
+}
+
+AliL3HoughTransformerRow::~AliL3HoughTransformerRow()
+{
+  DeleteHistograms();
+#ifdef do_mc
+  if(fTrackID)
+    {
+      for(Int_t i=0; i<GetNEtaSegments(); i++)
+       {
+         if(!fTrackID[i]) continue;
+         delete fTrackID[i];
+       }
+      delete [] fTrackID;
+      fTrackID = 0;
+    }
+#endif
+
+  if(fRowCount)
+    {
+      for(Int_t i=0; i<GetNEtaSegments(); i++)
+       {
+         if(!fRowCount[i]) continue;
+         delete [] fRowCount[i];
+       }
+      delete [] fRowCount;
+      fRowCount = 0;
+    }
+  if(fGapCount)
+    {
+      for(Int_t i=0; i<GetNEtaSegments(); i++)
+       {
+         if(!fGapCount[i]) continue;
+         delete [] fGapCount[i];
+       }
+      delete [] fGapCount;
+      fGapCount = 0;
+    }
+  if(fCurrentRowCount)
+    {
+      for(Int_t i=0; i<GetNEtaSegments(); i++)
+       {
+         if(fCurrentRowCount[i])
+           delete [] fCurrentRowCount[i];
+       }
+      delete [] fCurrentRowCount;
+      fCurrentRowCount = 0;
+    }
+}
+
+void AliL3HoughTransformerRow::DeleteHistograms()
+{
+  delete[] fLUT2sinphi0up;
+  delete[] fLUT2cosphi0up;
+  delete[] fLUT2sinphi0low;
+  delete[] fLUT2cosphi0low;
+
+  fLUT2sinphi0up=0;
+  fLUT2cosphi0up=0;
+  fLUT2sinphi0low=0;
+  fLUT2cosphi0low=0;
+
+  delete[] fLUTforwardZ;
+  delete[] fLUTforwardZ2;
+  delete[] fLUTbackwardZ;
+  delete[] fLUTbackwardZ2;
+
+  fLUTforwardZ=0;
+  fLUTforwardZ2=0;
+  fLUTbackwardZ=0;
+  fLUTbackwardZ2=0;
+
+  if(!fParamSpace)
+    return;
+  for(Int_t i=0; i<GetNEtaSegments(); i++)
+    {
+      if(!fParamSpace[i]) continue;
+      delete fParamSpace[i];
+    }
+  delete [] fParamSpace;
+}
+
+void AliL3HoughTransformerRow::CreateHistograms(Int_t nxbin,Float_t pt_min,
+                                               Int_t nybin,Float_t phimin,Float_t phimax)
+{
+  //Create the histograms (parameter space).
+  //These are 2D histograms, span by kappa (curvature of track) 
+  //and phi0 (emission angle with x-axis).
+  //The arguments give the range and binning; 
+  //nxbin = #bins in kappa
+  //nybin = #bins in phi0
+  //pt_min = mimium Pt of track (corresponding to maximum kappa)
+  //phi_min = mimimum phi0 
+  //phi_max = maximum phi0 
+    
+  Double_t x = AliL3Transform::GetBFact()*AliL3Transform::GetBField()/pt_min;
+  //Double_t torad = AliL3Transform::Pi()/180;
+  CreateHistograms(nxbin,-1.*x,x,nybin,phimin/**torad*/,phimax/**torad*/);
+}
+
+void AliL3HoughTransformerRow::CreateHistograms(Int_t nxbin,Float_t xmin,Float_t xmax,
+                                               Int_t nybin,Float_t ymin,Float_t ymax)
+{
+  
+  fParamSpace = new AliL3Histogram*[GetNEtaSegments()];
+  
+  Char_t histname[256];
+  for(Int_t i=0; i<GetNEtaSegments(); i++)
+    {
+      sprintf(histname,"paramspace_%d",i);
+      fParamSpace[i] = new AliL3Histogram(histname,"",nxbin,xmin,xmax,nybin,ymin,ymax);
+    }
+  
+#ifdef do_mc
+  if(fDoMC)
+    {
+      AliL3Histogram *hist = fParamSpace[0];
+      Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
+      if(!fTrackID)
+       {
+         cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(TrackIndex)<<" bytes to fTrackID"<<endl;
+         fTrackID = new TrackIndex*[GetNEtaSegments()];
+         for(Int_t i=0; i<GetNEtaSegments(); i++)
+           fTrackID[i] = new TrackIndex[ncells];
+       }
+    }
+#endif
+  AliL3Histogram *hist = fParamSpace[0];
+  Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
+  if(!fRowCount)
+    {
+      cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fRowCount"<<endl;
+      fRowCount = new UChar_t*[GetNEtaSegments()];
+      for(Int_t i=0; i<GetNEtaSegments(); i++)
+       fRowCount[i] = new UChar_t[ncells];
+    }
+  if(!fGapCount)
+    {
+      cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fGapCount"<<endl;
+      fGapCount = new UChar_t*[GetNEtaSegments()];
+      for(Int_t i=0; i<GetNEtaSegments(); i++)
+       fGapCount[i] = new UChar_t[ncells];
+    }
+  if(!fCurrentRowCount)
+    {
+      cout<<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fCurrentRowCount"<<endl;
+      fCurrentRowCount = new UChar_t*[GetNEtaSegments()];
+      for(Int_t i=0; i<GetNEtaSegments(); i++)
+       fCurrentRowCount[i] = new UChar_t[ncells];
+    }
+
+  //create lookup table for sin and cos
+  Int_t nbins = hist->GetNbinsY()+2;
+  fLUT2sinphi0up=new Float_t[nbins];
+  fLUT2cosphi0up=new Float_t[nbins];
+  fLUT2sinphi0low=new Float_t[nbins];
+  fLUT2cosphi0low=new Float_t[nbins];
+  Float_t hist_bin=hist->GetBinWidthY()/2.0;
+  for(Int_t i=hist->GetFirstYbin(); i<=hist->GetLastYbin(); i++){
+    Float_t phi0=hist->GetBinCenterY(i);
+    fLUT2sinphi0low[i]=2.*sin(phi0+hist_bin);
+    fLUT2cosphi0low[i]=2.*cos(phi0+hist_bin);
+    fLUT2sinphi0up[i]=2.*sin(phi0-hist_bin);
+    fLUT2cosphi0up[i]=2.*cos(phi0-hist_bin);
+  }  
+
+  //create lookup table for z of the digits
+  Int_t ntimebins = AliL3Transform::GetNTimeBins();
+  fLUTforwardZ = new Float_t[ntimebins];
+  fLUTforwardZ2 = new Float_t[ntimebins];
+  fLUTbackwardZ = new Float_t[ntimebins];
+  fLUTbackwardZ2 = new Float_t[ntimebins];
+  for(Int_t i=0; i<ntimebins; i++){
+    Float_t z;
+    z=AliL3Transform::GetZFast(0,i,GetZVertex());
+    fLUTforwardZ[i]=z;
+    fLUTforwardZ2[i]=z*z;
+    z=AliL3Transform::GetZFast(18,i,GetZVertex());
+    fLUTbackwardZ[i]=z;
+    fLUTbackwardZ2[i]=z*z;
+  }
+
+}
+
+void AliL3HoughTransformerRow::Reset()
+{
+  //Reset all the histograms
+
+  if(!fParamSpace)
+    {
+      LOG(AliL3Log::kWarning,"AliL3HoughTransformer::Reset","Histograms")
+       <<"No histograms to reset"<<ENDLOG;
+      return;
+    }
+  
+  for(Int_t i=0; i<GetNEtaSegments(); i++)
+    fParamSpace[i]->Reset();
+  
+#ifdef do_mc
+  if(fDoMC)
+    {
+      AliL3Histogram *hist = fParamSpace[0];
+      Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
+      for(Int_t i=0; i<GetNEtaSegments(); i++)
+       memset(fTrackID[i],0,ncells*sizeof(TrackIndex));
+    }
+#endif
+  AliL3Histogram *hist = fParamSpace[0];
+  Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
+  for(Int_t i=0; i<GetNEtaSegments(); i++)
+    {
+      memset(fRowCount[i],0,ncells*sizeof(UChar_t));
+      memset(fGapCount[i],1,ncells*sizeof(UChar_t));
+      memset(fCurrentRowCount[i],160,ncells*sizeof(UChar_t));
+    }
+}
+
+Int_t AliL3HoughTransformerRow::GetEtaIndex(Double_t eta)
+{
+  //Return the histogram index of the corresponding eta. 
+
+  Double_t etaslice = (GetEtaMax() - GetEtaMin())/GetNEtaSegments();
+  Double_t index = (eta-GetEtaMin())/etaslice;
+  return (Int_t)index;
+}
+
+inline AliL3Histogram *AliL3HoughTransformerRow::GetHistogram(Int_t eta_index)
+{
+  if(!fParamSpace || eta_index >= GetNEtaSegments() || eta_index < 0)
+    return 0;
+  if(!fParamSpace[eta_index])
+    return 0;
+  return fParamSpace[eta_index];
+}
+
+Double_t AliL3HoughTransformerRow::GetEta(Int_t eta_index,Int_t slice)
+{
+  Double_t eta_slice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments();
+  Double_t eta=0;
+  eta=(Double_t)((eta_index+0.5)*eta_slice);
+  return eta;
+}
+
+struct EtaRow {
+  UChar_t start_pad;
+  UChar_t end_pad;
+  Bool_t found;
+  Float_t start_y;
+#ifdef do_mc
+  Int_t mc_labels[MaxTrack];
+#endif
+};
+
+void AliL3HoughTransformerRow::TransformCircle()
+{
+  Int_t n_eta_segments = GetNEtaSegments();
+  Double_t eta_min = GetEtaMin();
+  Double_t eta_slice = (GetEtaMax() - eta_min)/n_eta_segments;
+
+  Int_t lower_threshold = GetLowerThreshold();
+
+  //Assumes that all the etaslice histos are the same!
+  AliL3Histogram *h = fParamSpace[0];
+  Int_t first_bin = h->GetFirstYbin();
+  Int_t last_bin = h->GetLastYbin();
+  Float_t x_min = h->GetXmin();
+  Float_t x_max = h->GetXmax();
+  Float_t x_bin = (x_max-x_min)/h->GetNbinsX();
+  Int_t first_binx = h->GetFirstXbin()-1;
+  Int_t last_binx = h->GetLastXbin()+1;
+  Int_t nbinx = h->GetNbinsX()+2;
+
+  UChar_t last_pad;
+  Int_t last_eta_index;
+  EtaRow *eta_clust = new EtaRow[n_eta_segments];
+
+  AliL3DigitRowData *tempPt = GetDataPointer();
+  if(!tempPt)
+    {
+      LOG(AliL3Log::kError,"AliL3HoughTransformer::TransformCircle","Data")
+       <<"No input data "<<ENDLOG;
+      return;
+    }
+
+  Int_t ipatch = GetPatch();
+  Double_t pad_pitch = AliL3Transform::GetPadPitchWidth(ipatch);
+  Int_t islice = GetSlice();
+  Float_t *fLUTz;
+  Float_t *fLUTz2;
+  if(islice < 18) {
+    fLUTz = fLUTforwardZ;
+    fLUTz2 = fLUTforwardZ2;
+  }
+  else {
+    fLUTz = fLUTbackwardZ;
+    fLUTz2 = fLUTbackwardZ2;
+  }
+
+  Int_t npads_in_patch = AliL3Transform::GetNPads(AliL3Transform::GetLastRow(ipatch))+2;
+  Bool_t **IsPrevRow = new Bool_t*[n_eta_segments];
+  Int_t nrows_in_patch = AliL3Transform::GetLastRow(ipatch)-AliL3Transform::GetFirstRow(ipatch)+2;
+  for(Int_t i=0; i<n_eta_segments; i++) {
+    IsPrevRow[i] = new Bool_t[npads_in_patch*nrows_in_patch];
+    memset(IsPrevRow[i],1,npads_in_patch*nrows_in_patch*sizeof(Bool_t));
+  }
+
+  //Loop over the padrows:
+  for(UChar_t i=AliL3Transform::GetFirstRow(ipatch); i<=AliL3Transform::GetLastRow(ipatch); i++)
+    {
+      Int_t npads = AliL3Transform::GetNPads((Int_t)i)-1;
+
+      last_pad = 255;
+      //Flush eta clusters array
+      memset(eta_clust,0,n_eta_segments*sizeof(EtaRow));  
+
+      Float_t x = AliL3Transform::Row2X((Int_t)i);
+      Float_t x2 = x*x;
+      Float_t y,r2;
+
+      Int_t lrow = (i-AliL3Transform::GetFirstRow(ipatch))*npads_in_patch;
+
+      //Get the data on this padrow:
+      AliL3DigitData *digPt = tempPt->fDigitData;
+      if((Int_t)i != (Int_t)tempPt->fRow)
+       {
+         cerr<<"AliL3HoughTransform::TransformCircle : Mismatching padrow numbering "<<(Int_t)i<<" "<<(Int_t)tempPt->fRow<<endl;
+         continue;
+       }
+      //      cout<<" Starting row "<<i<<endl;
+      //Loop over the data on this padrow:
+      for(UInt_t j=0; j<tempPt->fNDigit; j++)
+       {
+         UShort_t charge = digPt[j].fCharge;
+         if((Int_t)charge <= lower_threshold)
+           continue;
+         UChar_t pad = digPt[j].fPad;
+         UShort_t time = digPt[j].fTime;
+         if(pad != last_pad)
+           {
+             y = (pad-0.5*npads)*pad_pitch;
+             Float_t y2 = y*y;
+             r2 = x2 + y2;
+             last_eta_index = -1;
+           }
+
+         //Transform data to local cartesian coordinates:
+         Float_t z = fLUTz[(Int_t)time];
+         Float_t z2 = fLUTz2[(Int_t)time];
+         //Calculate the eta:
+         Double_t r = sqrt(r2+z2);
+         Double_t eta = 0.5 * log((r+z)/(r-z));
+         //Get the corresponding index, which determines which histogram to fill:
+         Int_t eta_index = (Int_t)((eta-eta_min)/eta_slice);
+
+#ifndef do_mc
+         if(eta_index == last_eta_index) continue;
+#endif
+         //      cout<<" Digit at patch "<<ipatch<<" row "<<i<<" pad "<<(Int_t)pad<<" time "<<time<<" etaslice "<<eta_index<<endl;
+         
+         if(eta_index < 0 || eta_index >= n_eta_segments)
+           continue;
+
+         if(!eta_clust[eta_index].found)
+           {
+             eta_clust[eta_index].start_pad = pad;
+             eta_clust[eta_index].end_pad = pad;
+             eta_clust[eta_index].start_y = y - pad_pitch/2.0;
+             eta_clust[eta_index].found = 1;
+#ifdef do_mc
+             if(fDoMC)
+               {
+                 for(Int_t t=0; t<3; t++)
+                   {
+                     Int_t label = digPt[j].fTrackID[t];
+                     if(label < 0) break;
+                     UInt_t c;
+                     for(c=0; c<(MaxTrack-1); c++)
+                       if(eta_clust[eta_index].mc_labels[c] == label || eta_clust[eta_index].mc_labels[c] == 0)
+                         break;
+                     //                      if(c == MaxTrack-1) cerr<<"AliL3HoughTransformer::TransformCircle : Cluster array reached maximum!! "<<c<<endl;
+                     eta_clust[eta_index].mc_labels[c] = label;
+                   }
+               }
+#endif
+             continue;
+           }
+         else
+           {
+             if(pad <= (eta_clust[eta_index].end_pad + 1))
+               {
+                 eta_clust[eta_index].end_pad = pad;
+#ifdef do_mc
+                 if(fDoMC)
+                   {
+                     for(Int_t t=0; t<3; t++)
+                       {
+                         Int_t label = digPt[j].fTrackID[t];
+                         if(label < 0) break;
+                         UInt_t c;
+                         for(c=0; c<(MaxTrack-1); c++)
+                           if(eta_clust[eta_index].mc_labels[c] == label || eta_clust[eta_index].mc_labels[c] == 0)
+                             break;
+                         //                      if(c == MaxTrack-1) cerr<<"AliL3HoughTransformer::TransformCircle : Cluster array reached maximum!! "<<c<<endl;
+                         eta_clust[eta_index].mc_labels[c] = label;
+                       }
+                   }
+#endif
+               }
+             else
+               {
+
+                 Bool_t fill_cluster = kFALSE;
+                 Bool_t *IsPrevRow2 = &IsPrevRow[eta_index][lrow+eta_clust[eta_index].start_pad];
+                 for(Int_t ipad = 0; ipad <= (eta_clust[eta_index].end_pad - eta_clust[eta_index].start_pad + 2); ipad++)
+                   {
+                     if(*IsPrevRow2)
+                       {
+                         fill_cluster = kTRUE;
+                         break;
+                       }
+                     IsPrevRow2++;
+                   }
+
+                 //              if(eta_index == 51)
+                 //                cout<<" Cluster to fill:"<<" row "<<(Int_t)i<<" pads from "<<(Int_t)eta_clust[eta_index].start_pad<<"("<<eta_clust[eta_index].start_y<<") to "<<(Int_t)eta_clust[eta_index].end_pad<<"("<<(eta_clust[eta_index].end_pad-0.5*(npads-1))*pad_pitch<<") fill "<<fill_cluster<<endl;
+
+                 if(fill_cluster) {
+                 //              cout<<" Cluster found at etaslice "<<eta_index<<" from pad "<<(Int_t)eta_clust[eta_index].start_pad<<" to pad "<<(Int_t)eta_clust[eta_index].end_pad<<endl;
+                 //              Bool_t fill_next_row = kFALSE;
+
+                 UChar_t *nrows = fRowCount[eta_index];
+                 UChar_t *ngaps = fGapCount[eta_index];
+                 UChar_t *currentrow = fCurrentRowCount[eta_index];
+
+                 //Do the transformation:
+                 Float_t start_y = eta_clust[eta_index].start_y;
+                 if(eta_clust[eta_index].start_pad == 0)
+                   start_y -= 3.0*pad_pitch;
+                 Float_t R1 = x2 + start_y*start_y;
+                 Float_t x_over_R1 = x/R1;
+                 Float_t start_y_over_R1 = start_y/R1;
+                 Float_t end_y = (eta_clust[eta_index].end_pad-0.5*(npads-1))*pad_pitch;
+                 if(eta_clust[eta_index].end_pad == npads)
+                   end_y += 3.0*pad_pitch;
+                 Float_t R2 = x2 + end_y*end_y; 
+                 Float_t x_over_R2 = x/R2;
+                 Float_t end_y_over_R2 = end_y/R2;
+
+                 //Fill the histogram along the phirange
+                 for(Int_t b=first_bin; b<=last_bin; b++)
+                   {
+                     Float_t kappa1 = start_y_over_R1*fLUT2cosphi0low[b]-x_over_R1*fLUT2sinphi0low[b];
+                     if(kappa1>x_max) continue;
+                     Float_t kappa2 = end_y_over_R2*fLUT2cosphi0up[b]-x_over_R2*fLUT2sinphi0up[b];
+                     if(kappa2<x_min) break;
+
+                     Int_t binx1 = 1 + (Int_t)((kappa1-x_min)/x_bin);
+                     if(binx1<first_binx) binx1 = first_binx;
+                     Int_t binx2 = 1 + (Int_t)((kappa2-x_min)/x_bin);
+                     if(binx2>last_binx) binx2 = last_binx;
+                     //                      if(eta_index == 51)
+                     //                        cout<<"     phi bin "<<b<<" kappa1 "<<kappa1<<" kappa2 "<<kappa2<<" from bin "<<binx1<<" to "<<binx2<<endl;
+                     Int_t temp_bin = b*nbinx + binx1;
+                     UChar_t *nrows2 = nrows + temp_bin;
+                     UChar_t *ngaps2 = ngaps + temp_bin;
+                     UChar_t *currentrow2 = currentrow + temp_bin;
+                     for(Int_t bin=binx1;bin<=binx2;bin++)
+                       {
+                         //Assumes threshold about 32
+                         if(*ngaps2 < 3) {
+                           if(i != *currentrow2)
+                             {
+                               (*nrows2)++;
+                               if(i > (*currentrow2 + 1))
+                                 (*ngaps2)++;
+                               *currentrow2=i;
+                             }
+#ifdef do_mc
+                           if(fDoMC)
+                             {
+                               for(UInt_t t=0;t<(MaxTrack-1); t++)
+                                 {
+                                   Int_t label = eta_clust[eta_index].mc_labels[t];
+                                   if(label == 0) break;
+                                   UInt_t c;
+                                   Int_t temp_bin2 = b*nbinx + bin;
+                                   for(c=0; c<(MaxTrack-1); c++)
+                                     if(fTrackID[eta_index][temp_bin2].fLabel[c] == label || fTrackID[eta_index][temp_bin2].fNHits[c] == 0)
+                                       break;
+                                   //                              if(c == MaxTrack-1) cerr<<"AliL3HoughTransformer::TransformCircle : Array reached maximum!! "<<c<<endl;
+                                   fTrackID[eta_index][temp_bin2].fLabel[c] = label;
+                                   if(fTrackID[eta_index][temp_bin2].fCurrentRow[c] != i) {
+                                     fTrackID[eta_index][temp_bin2].fNHits[c]++;
+                                     fTrackID[eta_index][temp_bin2].fCurrentRow[c] = i;
+                                   }
+                                 }
+                             }
+#endif
+                         }
+                         nrows2++;
+                         ngaps2++;
+                         currentrow2++;
+                       }
+                   }
+                 //End of the transformation
+
+                 /*              if(!fill_next_row) {
+                   Bool_t *IsCurrentRow = &IsPrevRow[eta_index][lrow+npads_in_patch];
+                   memset(&IsCurrentRow[eta_clust[eta_index].start_pad + 1],0,eta_clust[eta_index].end_pad - eta_clust[eta_index].start_pad + 1);
+                   } */
+                 }
+                 else {
+                   Bool_t *IsCurrentRow = &IsPrevRow[eta_index][lrow+npads_in_patch+eta_clust[eta_index].start_pad+1];
+                   memset(IsCurrentRow,0,eta_clust[eta_index].end_pad - eta_clust[eta_index].start_pad + 1);
+                 }
+
+                 Bool_t *IsCurrentRow = &IsPrevRow[eta_index][lrow+npads_in_patch+eta_clust[eta_index].end_pad+2];
+                 memset(IsCurrentRow,0,pad - eta_clust[eta_index].end_pad - 1);
+
+                 eta_clust[eta_index].start_pad = pad;
+                 eta_clust[eta_index].end_pad = pad;
+                 eta_clust[eta_index].start_y = y - pad_pitch/2.0;
+
+#ifdef do_mc
+                 if(fDoMC)
+                   {
+                     memset(eta_clust[eta_index].mc_labels,0,MaxTrack);
+                     for(Int_t t=0; t<3; t++)
+                       {
+                         Int_t label = digPt[j].fTrackID[t];
+                         if(label < 0) break;
+                         UInt_t c;
+                         for(c=0; c<(MaxTrack-1); c++)
+                           if(eta_clust[eta_index].mc_labels[c] == label || eta_clust[eta_index].mc_labels[c] == 0)
+                             break;
+                         //                      if(c == MaxTrack-1) cerr<<"AliL3HoughTransformer::TransformCircle : Cluster array reached maximum!! "<<c<<endl;
+                         eta_clust[eta_index].mc_labels[c] = label;
+                       }
+                   }
+#endif
+               }
+           }
+         last_pad = pad;
+         last_eta_index = eta_index;
+       }
+      //Write remaining clusters
+      for(Int_t eta_index = 0;eta_index < n_eta_segments;eta_index++)
+       {
+         //Check for empty row
+         if((eta_clust[eta_index].start_pad == 0) && (eta_clust[eta_index].end_pad == 0)) continue;
+
+         UChar_t *nrows = fRowCount[eta_index];
+         UChar_t *ngaps = fGapCount[eta_index];
+         UChar_t *currentrow = fCurrentRowCount[eta_index];
+
+         //Do the transformation:
+         Float_t start_y = eta_clust[eta_index].start_y;
+         if(eta_clust[eta_index].start_pad == 0)
+           start_y -= 3.0*pad_pitch;
+         Float_t R1 = x2 + start_y*start_y; 
+         Float_t x_over_R1 = x/R1;
+         Float_t start_y_over_R1 = start_y/R1;
+         Float_t end_y = (eta_clust[eta_index].end_pad-0.5*(npads-1))*pad_pitch;
+         if(eta_clust[eta_index].end_pad == npads)
+           end_y += 3.0*pad_pitch;
+         Float_t R2 = x2 + end_y*end_y; 
+         Float_t x_over_R2 = x/R2;
+         Float_t end_y_over_R2 = end_y/R2;
+
+         //Fill the histogram along the phirange
+         for(Int_t b=first_bin; b<=last_bin; b++)
+           {
+             Float_t kappa1 = start_y_over_R1*fLUT2cosphi0low[b]-x_over_R1*fLUT2sinphi0low[b];
+             if(kappa1>x_max) continue;
+             Float_t kappa2 = end_y_over_R2*fLUT2cosphi0up[b]-x_over_R2*fLUT2sinphi0up[b];
+             if(kappa2<x_min) break;
+
+             Int_t binx1 = 1 + (Int_t)((kappa1-x_min)/x_bin);
+             if(binx1<first_binx) binx1 = first_binx;
+             Int_t binx2 = 1 + (Int_t)((kappa2-x_min)/x_bin);
+             if(binx2>last_binx) binx2 = last_binx;
+
+             Int_t temp_bin = b*nbinx + binx1;
+             UChar_t *nrows2 = nrows + temp_bin;
+             UChar_t *ngaps2 = ngaps + temp_bin;
+             UChar_t *currentrow2 = currentrow + temp_bin;
+             for(Int_t bin=binx1;bin<=binx2;bin++)
+               {
+                 //Assumes threshold about 32
+                 if(*ngaps2 < 3) {
+                   if(i != *currentrow2)
+                     {
+                       (*nrows2)++;
+                       if(i > (*currentrow2 + 1))
+                         (*ngaps2)++;
+                       *currentrow2=i;
+                     }
+#ifdef do_mc
+                   if(fDoMC)
+                     {
+                       for(UInt_t t=0;t<(MaxTrack-1); t++)
+                         {
+                           Int_t label = eta_clust[eta_index].mc_labels[t];
+                           if(label == 0) break;
+                           UInt_t c;
+                           Int_t temp_bin2 = b*nbinx + bin;
+                           for(c=0; c<(MaxTrack-1); c++)
+                             if(fTrackID[eta_index][temp_bin2].fLabel[c] == label || fTrackID[eta_index][temp_bin2].fNHits[c] == 0)
+                               break;
+                           //                      if(c == MaxTrack-1) cout<<"AliL3HoughTransformer::TransformCircle : Array reached maximum!! "<<c<<endl;
+                           fTrackID[eta_index][temp_bin2].fLabel[c] = label;
+                           if(fTrackID[eta_index][temp_bin2].fCurrentRow[c] != i) {
+                             fTrackID[eta_index][temp_bin2].fNHits[c]++;
+                             fTrackID[eta_index][temp_bin2].fCurrentRow[c] = i;
+                           }
+                         }
+                     }
+#endif
+                 }
+                 nrows2++;
+                 ngaps2++;
+                 currentrow2++;
+               }
+           }
+         //End of the transformation
+
+         if(eta_clust[eta_index].end_pad < npads)
+           {
+             Bool_t *IsCurrentRow = &IsPrevRow[eta_index][lrow+npads_in_patch+eta_clust[eta_index].end_pad+2];
+             memset(IsCurrentRow,0,npads - eta_clust[eta_index].end_pad);
+           }
+
+       }
+
+      //Move the data pointer to the next padrow:
+      AliL3MemHandler::UpdateRowPointer(tempPt);
+    }
+
+
+  for(Int_t i=0; i<n_eta_segments; i++)
+    {
+      delete [] IsPrevRow[i];
+    }
+  delete [] IsPrevRow;
+
+  delete [] eta_clust;
+}
+
+Int_t AliL3HoughTransformerRow::GetTrackID(Int_t eta_index,Double_t kappa,Double_t psi)
+{
+  if(!fDoMC)
+    {
+      cerr<<"AliL3HoughTransformer::GetTrackID : Flag switched off"<<endl;
+      return -1;
+    }
+  
+#ifdef do_mc
+  if(eta_index < 0 || eta_index > GetNEtaSegments())
+    {
+      cerr<<"AliL3HoughTransformer::GetTrackID : Wrong etaindex "<<eta_index<<endl;
+      return -1;
+    }
+  AliL3Histogram *hist = fParamSpace[eta_index];
+  Int_t bin = hist->FindBin(kappa,psi);
+  Int_t label=-1;
+  Int_t max=0;
+  for(UInt_t i=0; i<(MaxTrack-1); i++)
+    {
+      Int_t nhits=fTrackID[eta_index][bin].fNHits[i];
+      if(nhits == 0) break;
+      if(nhits > max)
+       {
+         max = nhits;
+         label = fTrackID[eta_index][bin].fLabel[i];
+       }
+    }
+  Int_t label2=-1;
+  Int_t max2=0;
+  for(UInt_t i=0; i<(MaxTrack-1); i++)
+    {
+      Int_t nhits=fTrackID[eta_index][bin].fNHits[i];
+      if(nhits == 0) break;
+      if(nhits > max2)
+       {
+         if(fTrackID[eta_index][bin].fLabel[i]!=label) {
+           max2 = nhits;
+           label2 = fTrackID[eta_index][bin].fLabel[i];
+         }
+       }
+    }
+  //  cout<<" TrackID label "<<label<<" max "<<max<<" label2 "<<label2<<" max2 "<<max2<<" "<<(Float_t)max2/(Float_t)max<<" "<<fTrackID[eta_index][bin].fLabel[MaxTrack-1]<<" "<<(Int_t)fTrackID[eta_index][bin].fNHits[MaxTrack-1]<<endl;
+  return label;
+#endif
+  cout<<"AliL3HoughTransformer::GetTrackID : Compile with do_mc flag!"<<endl;
+  return -1;
+}
+
+UChar_t *AliL3HoughTransformerRow::GetRowCount(Int_t eta_index)
+{
+  return fRowCount[eta_index];
+}
+
+UChar_t *AliL3HoughTransformerRow::GetGapCount(Int_t eta_index)
+{
+  return fGapCount[eta_index];
+}
similarity index 54%
rename from HLT/hough/AliL3HoughTransformerGap.h
rename to HLT/hough/AliL3HoughTransformerRow.h
index b975f78..159bbaf 100644 (file)
@@ -1,36 +1,26 @@
 // @(#) $Id$
 
-#ifndef ALIL3_HOUGHTRANSFORMERGAP
-#define ALIL3_HOUGHTRANSFORMERGAP
+#ifndef ALIL3_HOUGHTRANSFORMERROW
+#define ALIL3_HOUGHTRANSFORMERROW
 
 #include "AliL3RootTypes.h"
 #include "AliL3HoughBaseTransformer.h"
 
 class AliL3Histogram;
 
-#define COUNTINGROWS
-//#define SUBSLICES
-
-class AliL3HoughTransformerGap : public AliL3HoughBaseTransformer {
+class AliL3HoughTransformerRow : public AliL3HoughBaseTransformer {
   
  private:
   
   AliL3Histogram **fParamSpace; //!
 
 #ifdef do_mc
-  TrackIndex **fTrackID; //!
+  static TrackIndex **fTrackID; //!
 #endif
   Bool_t fDoMC;
 
-  Float_t fZVertex;
-
   void DeleteHistograms();
 
-#ifdef SUBSLICES
-  Int_t **fEtaSpace; //!
-  Int_t fEtaSpaceSize;
-#endif
-#ifdef COUNTINGROWS
   static UChar_t **fRowCount; //!
   static UChar_t **fGapCount; //!
   static UChar_t **fCurrentRowCount; //!
@@ -39,12 +29,16 @@ class AliL3HoughTransformerGap : public AliL3HoughBaseTransformer {
   Float_t *fLUT2cosphi0up; //!
   Float_t *fLUT2sinphi0low; //!   
   Float_t *fLUT2cosphi0low; //!
-#endif
+
+  Float_t *fLUTforwardZ; //!
+  Float_t *fLUTforwardZ2; //!
+  Float_t *fLUTbackwardZ; //!
+  Float_t *fLUTbackwardZ2; //!
 
  public:
-  AliL3HoughTransformerGap(); 
-  AliL3HoughTransformerGap(Int_t slice,Int_t patch,Int_t n_eta_segments,Bool_t DoEtaOverlap=kFALSE,Bool_t DoMC=kFALSE,Float_t zvertex=0.0);
-  virtual ~AliL3HoughTransformerGap();
+  AliL3HoughTransformerRow(); 
+  AliL3HoughTransformerRow(Int_t slice,Int_t patch,Int_t n_eta_segments,Bool_t DoMC=kFALSE,Float_t zvertex=0.0);
+  virtual ~AliL3HoughTransformerRow();
 
   //void CreateHistograms(Float_t ptmin,Float_t ptmax,Float_t ptres,Int_t nybin,Float_t psi);
   void CreateHistograms(Int_t nxbin,Float_t ptmin,Int_t nybin,Float_t phimin,Float_t phimax);
@@ -56,22 +50,11 @@ class AliL3HoughTransformerGap : public AliL3HoughBaseTransformer {
   Int_t GetEtaIndex(Double_t eta);
   AliL3Histogram *GetHistogram(Int_t eta_index);
   Double_t GetEta(Int_t eta_index,Int_t slice);
-
-#ifdef SUBSLICES
-  Int_t GetEtaSubIndex(Double_t eta,Int_t eta_index);
-  Double_t GetMaxSubEta(Float_t kappa,Float_t phi0,Int_t eta_index,Int_t slice);
-#endif
-
-  //Int_t GetTrackID(Int_t eta_index,Double_t kappa,Double_t psi);
-
+  Int_t GetTrackID(Int_t eta_index,Double_t kappa,Double_t psi);
   UChar_t *GetRowCount(Int_t eta_index);
   UChar_t *GetGapCount(Int_t eta_index);
 
-#ifdef COUNTINGROWS
-  //void CalculateNRows();
-#endif
-
-  ClassDef(AliL3HoughTransformerGap,1) //Hough transformation class
+  ClassDef(AliL3HoughTransformerRow,1) //TPC Rows Hough transformation class
 
 };
 
index 57b617d..121ef86 100644 (file)
@@ -13,7 +13,7 @@
 #include "AliL3HoughTransformerVhdl.h"
 #include "AliL3FFloat.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 041d479..6a70b28 100644 (file)
@@ -18,7 +18,7 @@ SRCS = AliL3HoughTransformer.cxx AliL3HoughClusterTransformer.cxx \
        AliL3HoughEval.cxx AliL3HoughMerger.cxx AliL3HoughBaseTransformer.cxx \
        AliL3HoughIntMerger.cxx AliL3HoughGlobalMerger.cxx AliL3HoughTransformerVhdl.cxx \
        AliL3Histogram.cxx AliL3Histogram1D.cxx AliL3HoughMaxFinder.cxx AliL3Hough.cxx \
-       AliL3HoughTransformerLUT.cxx AliL3HoughTransformerGap.cxx AliL3HistogramAdaptive.cxx 
+       AliL3HoughTransformerLUT.cxx AliL3HoughTransformerRow.cxx AliL3HistogramAdaptive.cxx 
 
 ifeq ($(ARCH),Darwin)
 ## AliL3HoughTrack put into src as symbolic link
index 54a793c..6098508 100644 (file)
@@ -54,7 +54,6 @@ AliL3Kalman::AliL3Kalman(Char_t *datapath, Int_t *slice, Int_t min_clusters = 0)
   fWriteOut = kTRUE;
   fBenchmark = 0;
 
-  fMakeSeed = kFALSE;
 }  
 
 AliL3Kalman::~AliL3Kalman()
@@ -153,33 +152,10 @@ void AliL3Kalman::ProcessTracks()
 
   fKalmanTracks = new AliL3TrackArray();
 
-  // Make a tree to store state vector, covariance matrix and chisquare
+  // Make a ntuple to store state vector, covariance matrix and chisquare
   // Will eventually not need a TTree??
-  TTree *kalmanTree = new TTree("kalmanTree","kalmantracks");
+  TNtuple *kalmanTree = new TNtuple("kalmanTree","kalmantracks","x0:x1:x2:x3:x4:c0:c1:c2:c3:c4:c5:c6:c7:c8:c9:c10:c11:c12:c13:c14:chisq");
   Float_t meas[21];
-  Int_t lab = 123456789;
-  kalmanTree->Branch("x0",&meas[0],"x0/F");
-  kalmanTree->Branch("x1",&meas[1],"x1/F");
-  kalmanTree->Branch("x2",&meas[2],"x2/F");
-  kalmanTree->Branch("x3",&meas[3],"x3/F");
-  kalmanTree->Branch("x4",&meas[4],"x4/F");
-  kalmanTree->Branch("c0",&meas[5],"c0/F");
-  kalmanTree->Branch("c1",&meas[6],"c1/F");
-  kalmanTree->Branch("c2",&meas[7],"c2/F");
-  kalmanTree->Branch("c3",&meas[8],"c3/F");
-  kalmanTree->Branch("c4",&meas[9],"c4/F");
-  kalmanTree->Branch("c5",&meas[10],"c5/F");
-  kalmanTree->Branch("c6",&meas[11],"c6/F");
-  kalmanTree->Branch("c7",&meas[12],"c7/F");
-  kalmanTree->Branch("c8",&meas[13],"c8/F");
-  kalmanTree->Branch("c9",&meas[14],"c9/F");
-  kalmanTree->Branch("c10",&meas[15],"c10/F");
-  kalmanTree->Branch("c11",&meas[16],"c11/F");
-  kalmanTree->Branch("c12",&meas[17],"c12/F");
-  kalmanTree->Branch("c13",&meas[18],"c13/F");
-  kalmanTree->Branch("c14",&meas[19],"c14/F");
-  kalmanTree->Branch("chisq",&meas[20],"chisq/F");
-  kalmanTree->Branch("lab",&lab,"lab/I");
 
   // Go through the tracks from conformal or hough tracker
   for (Int_t iTrack = 0; iTrack < fTracks->GetNTracks(); iTrack++)
@@ -192,41 +168,20 @@ void AliL3Kalman::ProcessTracks()
       if (!track) continue;
       if (track->GetNumberOfPoints() < fMinPointsOnTrack) continue;    
 
-      UInt_t *hitnum = track->GetHitNumbers();
-      UInt_t id;
-      
-      id = hitnum[0];
-      Int_t slice = (id>>25) & 0x7f;
-      Int_t patch = (id>>22) & 0x7;    
-      UInt_t pos = id&0x3fffff;
-      AliL3SpacePointData *points = fClusters[slice][patch];
-      lab=points[pos].fTrackID[0];
-      //cout << lab << endl;
-
       AliL3KalmanTrack *kalmantrack = new AliL3KalmanTrack();
 
       Bool_t save = kTRUE;
-      if (fMakeSeed == kTRUE)
-       {
-         if (MakeSeed(kalmantrack, track) == 0)
-           {
-             save = kFALSE;
-             continue;
-           }
-       }
-      else 
+
+      if (InitKalmanTrack(kalmantrack, track) == 0)
        {
-         if (InitKalmanTrack(kalmantrack, track) == 0)
-           {
-             save = kFALSE;
-             continue;
-           }
+         save = kFALSE;
+         continue;
        }
+      
 
       if (Propagate(kalmantrack, track) == 0) 
        {
          save = kFALSE;
-         continue;
        }
       
       if (save) {
@@ -262,8 +217,8 @@ void AliL3Kalman::ProcessTracks()
        outtrack->Set(track);
 
        // Fill the ntuple with the state vector, covariance matrix and
-       // chisquare and track label
-       kalmanTree->Fill();
+       // chisquare
+       kalmanTree->Fill(meas);
       }
 
       delete track;
@@ -311,39 +266,7 @@ Int_t AliL3Kalman::InitKalmanTrack(AliL3KalmanTrack *kalmantrack, AliL3Track *tr
   UInt_t pos = id&0x3fffff;
   AliL3SpacePointData *points = fClusters[slice][patch];
 
-  return kalmantrack->Init(track, points, pos);
-}
-
-Int_t AliL3Kalman::MakeSeed(AliL3KalmanTrack *kalmantrack, AliL3Track *track)
-{  
-  // Makes a rough state vector and covariance matrix based on three
-  // space points of the loaded track 
-  // Will eventually be removed?? Will use track parameters from conformal
-  // tracker or HT as seeds for the kalman filter.
-  UInt_t *hitnum = track->GetHitNumbers();
-  UInt_t id1, id2, id3;
-
-  id1 = hitnum[0];
-  Int_t slice1 = (id1>>25) & 0x7f;
-  Int_t patch1 = (id1>>22) & 0x7;      
-  UInt_t pos1 = id1&0x3fffff;
-  AliL3SpacePointData *points1 = fClusters[slice1][patch1];
-  
-  id2 = hitnum[Int_t(track->GetNHits()/2)];
-  Int_t slice2 = (id2>>25) & 0x7f;
-  Int_t patch2 = (id2>>22) & 0x7;      
-  UInt_t pos2 = id2&0x3fffff;
-  AliL3SpacePointData *points2 = fClusters[slice2][patch2];
-
-  id3 = hitnum[track->GetNHits()-1];
-  Int_t slice3 = (id3>>25) & 0x7f;
-  Int_t patch3 = (id3>>22) & 0x7;      
-  UInt_t pos3 = id3&0x3fffff;
-  AliL3SpacePointData *points3 = fClusters[slice3][patch3];
-
-  return kalmantrack->MakeTrackSeed(points1,pos1,points2,pos2,points3,pos3);
-
+  return kalmantrack->Init(track, points, pos, slice);
 }
 
 Int_t AliL3Kalman::Propagate(AliL3KalmanTrack *kalmantrack, AliL3Track *track)
@@ -365,7 +288,7 @@ Int_t AliL3Kalman::Propagate(AliL3KalmanTrack *kalmantrack, AliL3Track *track)
 
       AliL3SpacePointData *points = fClusters[slice][patch];
       if (!points) continue;
-      if (kalmantrack->Propagate(points,pos) == 0) 
+      if (kalmantrack->Propagate(points,pos,slice) == 0) 
        {
          badpoint++;
          continue;
index ee4f71d..2d6f437 100644 (file)
@@ -31,21 +31,18 @@ class AliL3Kalman {
   Char_t fWriteOutPath[256];
   Bool_t fWriteOut;
   Int_t fEvent;
-  Bool_t fMakeSeed;
 
  public:
 
-  AliL3Kalman(Char_t *datapath, Int_t *slice=0, Int_t min_clusters=0);
+  AliL3Kalman(Char_t *datapath, Int_t *slice, Int_t min_clusters);
   virtual ~AliL3Kalman();
   void Init();
   void LoadTracks(Int_t event, Bool_t sp);
   void ProcessTracks();
   Int_t InitKalmanTrack(AliL3KalmanTrack *kalmantrack, AliL3Track *track);
-  Int_t MakeSeed(AliL3KalmanTrack *kalmantrack, AliL3Track *track);
   Int_t Propagate(AliL3KalmanTrack *kalmantrack, AliL3Track *track);
   Int_t Update(AliL3SpacePointData *points, UInt_t pos, AliL3KalmanTrack *kalmantrack);
   void WriteFiles(Char_t *path="data"){fWriteOut = kTRUE; sprintf(fWriteOutPath,"%s",path);}
-  void DoMakeSeed(){fMakeSeed = kTRUE;}
   Double_t GetCpuTime();
   AliL3TrackArray *GetTracks() {return fKalmanTracks;}
 };
index 458e634..bfea605 100644 (file)
@@ -3,12 +3,23 @@
 #include "AliL3SpacePointData.h"
 #include "AliL3Logging.h"
 #include "AliL3Transform.h"
+#include "AliL3StandardIncludes.h"
+
+// includes for offline comparison, will be removed
+#include "AliTPCtrack.h"
+// includes for offline comparison, will be removed
+
+#include "Riostream.h" 
+
+
 ClassImp(AliL3KalmanTrack)
 
 // Class for kalman tracks
 
 AliL3KalmanTrack::AliL3KalmanTrack()
 {
+  fX = 0;
+
   fMaxChi2 = 12;
   // Constructor
 }
@@ -18,30 +29,50 @@ AliL3KalmanTrack::~AliL3KalmanTrack()
   // Destructor
 }
 
-Int_t AliL3KalmanTrack::Init(AliL3Track *track, AliL3SpacePointData *points, UInt_t pos)
+Int_t AliL3KalmanTrack::Init(AliL3Track *track, AliL3SpacePointData *points, UInt_t pos, Int_t slice)
 {
-  fX = points[pos].fX; 
 
-  fP0 = points[pos].fY;
+  Float_t xyz[3];
+
+  xyz[0] = points[pos].fX;
+  xyz[1] = points[pos].fY;
+
+  // NB! I think boolean variable in the command under is true if single slice.
+  // Better do a test if it is signle slice, and if so set it true.??.
+  AliL3Transform::Global2Local(xyz,slice);
+  
+  fX = xyz[0];
+
+  fP0 = xyz[1];
   fP1 = points[pos].fZ; 
   fP3 = track->GetTgl(); //cout << fP3; 
   if (TMath::Abs(fP3) > 1.2) return 0; //From AliTPCtrackerMI
   fP4 = 0.5*(-track->GetCharge()*1./track->GetPt())*0.0029980*AliL3Transform::GetBField(); // Calculation of the curvature 
-  
   if (TMath::Abs(fP4) >= 0.0066) return 0; // From AliTPCtrackerMI
 
-  Float_t centerX = track->GetFirstPointX() - ((track->GetPt()/(0.0029980*AliL3Transform::GetBField())) * TMath::Cos(track->GetPsi() + track->GetCharge() * 0.5 * 3.14159265358979323846));
-  Float_t centerY = track->GetFirstPointY() - ((track->GetPt()/(0.0029980*AliL3Transform::GetBField())) * TMath::Sin(track->GetPsi() + track->GetCharge() * 0.5 * 3.14159265358979323846));
+
+  Float_t firstXY[2];
+  firstXY[0] = track->GetFirstPointX();
+  firstXY[1] = track->GetFirstPointY();
+
+  AliL3Transform::Global2Local(firstXY,slice);
+
+  //cout << firstXY[0] << endl;
+
+  //Float_t centerX = track->GetFirstPointX() - ((track->GetPt()/(0.0029980*AliL3Transform::GetBField())) * TMath::Cos(track->GetPsi() + track->GetCharge() * 0.5 * 3.14159265358979323846));
+  Float_t centerX = firstXY[0] - ((track->GetPt()/(0.0029980*AliL3Transform::GetBField())) * TMath::Cos(track->GetPsi() + track->GetCharge() * 0.5 * 3.14159265358979323846));
+  //Float_t centerY = track->GetFirstPointY() - ((track->GetPt()/(0.0029980*AliL3Transform::GetBField())) * TMath::Sin(track->GetPsi() + track->GetCharge() * 0.5 * 3.14159265358979323846));
+
   fP2 = fP4*centerX; // Curvature times center of curvature
-  
+
   if (TMath::Abs(fP4*fX - fP2) >= 0.9) // What's this??
     {
       return 0;
     }
   //cout << track->GetCharge() << endl; 
   //cout << "AliL3KalmanTrack::Init, " << fP0 << " " << fP1 << " " << fP2 << " " << fP3 << " " << fP4 << endl;
-  //cout << fP2 << " " << fP4 << endl;
-  Float_t num = 12;
+  //cout << fP4 << endl;
+  //Float_t num = 12;
 
   fC00 = points[pos].fSigmaY2; 
   fC10 = 0; fC11 = points[pos].fSigmaZ2;
@@ -60,85 +91,12 @@ Int_t AliL3KalmanTrack::Init(AliL3Track *track, AliL3SpacePointData *points, UIn
   return 1;
 }
 
-Int_t AliL3KalmanTrack::MakeTrackSeed(AliL3SpacePointData *points1, UInt_t pos1, AliL3SpacePointData *points2, UInt_t pos2, AliL3SpacePointData *points3, UInt_t pos3)
-{
-  // Make track seed based on three clusters 
-  fChisq = 0;
-  fX = points1[pos1].fX;
-
-  fP0 = points1[pos1].fY; 
-  fP1 = points1[pos1].fZ; 
-
-  Float_t X2 = points2[pos2].fX;
-  Float_t Y2 = points2[pos2].fY;
-  Float_t Z2 = points2[pos2].fZ;
-
-  Float_t X3 = points3[pos3].fX;
-  Float_t Y3 = points3[pos3].fY;
-  Float_t Z3 = points3[pos3].fZ;
-
-  Float_t ZZ = fP1 - ((fP1 - Z3)/(fX-X3))*(fX-X2);
-  if (TMath::Abs(ZZ - Z2) > 10) return 0; //What's this?? (fP1 - Z3)/(fX-X3)*(fX-X2) is an angle
-
-  // It may make no difference. Check on a big event??.
-  if ((X2-fX)*(0-Y2)-(0-X2)*(Y2-fP0) == 0) return 0; //Straight seed
-
-  // Initial approximation of the state vector
-  fP2 = f2(fX,fP0,X2,Y2,X3,Y3);
-  fP3 = f3(fX,fP0,X2,Y2,fP1,Z2);
-  fP4 = f4(fX,fP0,X2,Y2,X3,Y3);
-  if (TMath::Abs(fP3) > 1.2) return 0; //From AliTPCtrackerMI
-  if (TMath::Abs(fP4) >= 0.0066) return 0; // From AliTPCtrackerMI
-
-  //cout << "AliL3KalmanTrack::MakeTrackSeed, " << fP0 << " " << fP1 << " " << fP2 << " " << fP3 << " " << fP4 << endl;
-  //cout << fP2 << " " << fP4 << endl;
-
-  // Initial approximation of the covariance matrix
-  Float_t sy1=points1[pos1].fSigmaY2;
-  Float_t sz1=points1[pos1].fSigmaZ2;
-  Float_t sy2=points2[pos2].fSigmaY2;
-  Float_t sz2=points2[pos2].fSigmaZ2;
-  Float_t sy3=25000*fP4*fP4+0.1;
-  Float_t sy=0.1;
-  Float_t sz=0.1;
-  
-  Float_t f40=(f4(fX,fP0+sy,X2,Y2,X3,Y3)-fP4)/sy;
-  Float_t f42=(f4(fX,fP0,X2,Y2+sy,X3,Y3)-fP4)/sy;
-  Float_t f43=(f4(fX,fP0,X2,Y2,X3,Y3+sy)-fP4)/sy;
-  Float_t f20=(f2(fX,fP0+sy,X2,Y2,X3,Y3)-fP2)/sy;
-  Float_t f22=(f2(fX,fP0,X2,Y2+sy,X3,Y3)-fP2)/sy;
-  Float_t f23=(f2(fX,fP0,X2,Y2,X3,Y3+sy)-fP2)/sy;
-  Float_t f30=(f3(fX,fP0+sy,X2,Y2,fP1,Z2)-fP3)/sy;
-  Float_t f31=(f3(fX,fP0,X2,Y2,fP1+sz,Z2)-fP3)/sz;
-  Float_t f32=(f3(fX,fP0,X2,Y2+sy,fP1,Z2)-fP3)/sy;
-  Float_t f34=(f3(fX,fP0,X2,Y2,fP1,Z2+sz)-fP3)/sz;
-  
-  fC00 = sy1;
-  fC10 = 0;  
-  fC11 = sz1;
-  fC20 = f20*sy1;  
-  fC21 = 0;  
-  fC22 = f20*sy1*f20+f22*sy2*f22+f23*sy3*f23;
-  fC30 = f30*sy1;  
-  fC31 = f31*sz1;  
-  fC32 = f30*sy1*f20+f32*sy2*f22;
-  fC33 = f30*sy1*f30+f31*sz1*f31+f32*sy2*f32+f34*sz2*f34;
-  fC40 = f40*sy1; 
-  fC41 = 0; 
-  fC42 = f40*sy1*f20+f42*sy2*f22+f43*sy3*f23;
-  fC43 = f30*sy1*f40+f32*sy2*f42; 
-  fC44 = f40*sy1*f40+f42*sy2*f42+f43*sy3*f43;
-  //cout << "Errors: " << fC00 << " " << fC11 << " " << fC20 << " " << fC22 << " " << " " << fC30 << " " << fC31 << " " << fC32 << " " << fC33 << " " << fC40 << " " << fC42 << " " << fC43 << " " << fC44 << endl;
-  
-  return 1;
-}
-
-Int_t AliL3KalmanTrack::Propagate(AliL3SpacePointData *points, UInt_t pos)
+Int_t AliL3KalmanTrack::Propagate(AliL3SpacePointData *points, UInt_t pos, Int_t slice)
 {
   // Propagates track to the plane of the next found cluster
   Float_t Xold = fX; // X position for previous space point
-  Float_t Xnew = points[pos].fX; // X position of current space point
-  Float_t dx = Xnew - Xold;
+  //Float_t Xnew = points[pos].fX; // X position of current space point
+  //Float_t dx = Xnew - Xold; cout << Xnew << endl;
   Float_t Yold = fP0; // Y position of old point
   Float_t Zold = fP1; // Z position of old point
   Float_t par2old = fP2;  
@@ -160,6 +118,16 @@ Int_t AliL3KalmanTrack::Propagate(AliL3SpacePointData *points, UInt_t pos)
   Float_t oldfC43 = fC43;
   Float_t oldfC44 = fC44;
 
+  Float_t xyz[3];
+
+  xyz[0] = points[pos].fX;
+  xyz[1] = points[pos].fY;
+
+  AliL3Transform::Global2Local(xyz,slice);
+  
+  Float_t Xnew = xyz[0];
+  Float_t dx = Xnew - Xold; //cout << Xnew << endl;
+  
   if (TMath::Abs(fP4*Xnew - fP2) >= 0.9) // What's this??
     {
       //cout << "Propagation failed! Stiff track. " << pos << endl;
@@ -202,6 +170,7 @@ Int_t AliL3KalmanTrack::Propagate(AliL3SpacePointData *points, UInt_t pos)
   Float_t f02 = -dx*(2*RR + CC*(Cold/Rold + Cnew/Rnew))/(RR*RR);
   Float_t f04 = dx*(RR*XX + CC*(Cold*Xold/Rold + Cnew*Xnew/Rnew))/(RR*RR);
   Float_t CR = Cold*Rnew + Cnew*Rold;
+
   Float_t f12 = -dx*fP3*(2*CR + CC*(Cnew*(Cold/Rold)-Rold + Cold*(Cnew/Rnew)-Rnew))/(CR*CR);
   Float_t f13 = dx*CC/CR;
   Float_t f14 = dx*fP3*(CR*XX-CC*(Rold*Xnew-Cnew*Cold*Xold/Rold + Rnew*Xold-Cold*Cnew*Xnew/Rnew))/(CR*CR);
@@ -271,9 +240,10 @@ Int_t AliL3KalmanTrack::Propagate(AliL3SpacePointData *points, UInt_t pos)
   fP2 += Xnew*(fP4-CC);
   */
   // Update the track parameters with the measured values of the new point
-  Int_t value = UpdateTrack(points, pos);
 
-  /*if (value == 0)
+  Int_t value = UpdateTrack(points, pos, slice);
+
+  if (value == 0)
     {
       fP0 = Yold;
       fP1 = Zold;
@@ -297,17 +267,26 @@ Int_t AliL3KalmanTrack::Propagate(AliL3SpacePointData *points, UInt_t pos)
       fC44 = oldfC44;
 
       return value;
-      }*/
+      }
   
-  //else 
-return value;
+  else 
+    return value;
   //return UpdateTrack(points, pos);
 }
 
-Int_t AliL3KalmanTrack::UpdateTrack(AliL3SpacePointData *points, UInt_t pos)
+Int_t AliL3KalmanTrack::UpdateTrack(AliL3SpacePointData *points, UInt_t pos, Int_t slice)
 {
   // Update the track parameters with the measured values
-  fX = points[pos].fX; 
+  
+  Float_t xyz[3];
+
+  xyz[0] = points[pos].fX;
+  xyz[1] = points[pos].fY;
+  
+  AliL3Transform::Global2Local(xyz,slice);
+  
+  //fX = points[pos].fX; 
+  fX = xyz[0];
 
   // The errors from the measurement of the spacepoint
   Float_t sigmaY2 = points[pos].fSigmaY2;
@@ -337,9 +316,10 @@ Int_t AliL3KalmanTrack::UpdateTrack(AliL3SpacePointData *points, UInt_t pos)
   Float_t k41 = fC40*sigmaYZ + fC41*sigmaZ2;
 
   // Deviation between the predicted and measured values of y and z 
-  Float_t dy = points[pos].fY-fP0; //cout << "dy = " << dy;
+  //Float_t dy = points[pos].fY-fP0; //cout << "dy = " << dy;
+  Float_t dy = xyz[1] - fP0; //cout << "dy = " << dy;
   Float_t dz = points[pos].fZ-fP1; //cout << ", dz = " << dz << endl; 
-  //cout << "Measured " << points[pos].fY << " " << points[pos].fZ << endl;
+  //cout << "Measured " << xyz[2] << " " << points[pos].fZ << endl;
 
   // Prediction of fP2 and fP4
   Float_t cur = fP4 + k40*dy + k41*dz; 
@@ -392,20 +372,22 @@ Int_t AliL3KalmanTrack::UpdateTrack(AliL3SpacePointData *points, UInt_t pos)
   //cout << "AliL3KalmanTrack::Update, Chi2 = " << GetChisq() << endl;
   //cout << "AliKalmanTrack::Update, sigmaY2 = " << sigmaY2 << " sigmaZ2 = " << sigmaZ2 << " sigmaYZ = " << sigmaYZ << " dy = " << dy << " dz = " << dz << endl;
   // Calculate increase of chisquare
-  fChisq = GetChisq() + GetChisqIncrement(points,pos);
+  fChisq = GetChisq() + GetChisqIncrement(xyz[1],sigmaY2,points[pos].fZ,sigmaZ2);
   //cout << "fChisq = " << fChisq << endl;
 //(dy*sigmaY2*dy + 2*sigmaYZ*dy*dz + dz*sigmaZ2*dz) / (sigmaY2*sigmaZ2 - sigmaYZ*sigmaYZ);
   // Must at some point make an cut on chisq. Here?
-  //  if (fChisq > fMaxChi2) return 0;
+  //if (fChisq > fMaxChi2) return 0;
 
   return 1;
 } 
 
-Float_t AliL3KalmanTrack::GetChisqIncrement(AliL3SpacePointData *points, UInt_t pos)
+//Float_t AliL3KalmanTrack::GetChisqIncrement(AliL3SpacePointData *points, UInt_t pos)
+Float_t AliL3KalmanTrack::GetChisqIncrement(Float_t y, Float_t error_y, Float_t z, Float_t error_z)
 {
   // This function calculates a predicted chi2 increment.
   //-----------------------------------------------------------------
-  Float_t r00=points[pos].fSigmaY2, r01=0., r11=points[pos].fSigmaZ2;
+  //Float_t r00=points[pos].fSigmaY2, r01=0., r11=points[pos].fSigmaZ2;
+  Float_t r00=error_y, r01=0., r11=error_z;
   r00+=fC00; r01+=fC10; r11+=fC11;
 
   Double_t det=r00*r11 - r01*r01;
@@ -416,7 +398,7 @@ Float_t AliL3KalmanTrack::GetChisqIncrement(AliL3SpacePointData *points, UInt_t
     }*/
   Double_t tmp=r00; r00=r11; r11=tmp; r01=-r01;
   
-  Double_t dy=points[pos].fY - fP0, dz=points[pos].fZ - fP1;
+  Double_t dy=y - fP0, dz=z - fP1;
   //cout << "AliTPCtrack::GetPredictedChi2, r00 = " << r00 << " r11 = " << r11 << " ro1 = " << r01 << " dy = " << dy << " dz = " << dz << endl;
   return (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz)/det;
 }
@@ -498,3 +480,265 @@ void AliL3KalmanTrack::Set(AliL3KalmanTrack *track)
 
   SetNHits(tpt->GetNHits());
 }
+
+Int_t AliL3KalmanTrack::PropagateOfflineTrack(Double_t x, Double_t y, Double_t z, Double_t ey, Double_t ez)
+{
+  // Propagates track to the plane of the next found cluster
+  Float_t Xold = fX; // X position for previous space point
+  Float_t Xnew = x; // X position of current space point
+  Float_t dx = Xnew - Xold; //cout << dx << endl;
+  Float_t Yold = fP0; // Y position of old point
+  Float_t Zold = fP1; // Z position of old point
+  Float_t par2old = fP2;
+  Float_t par3old = fP3;
+  Float_t par4old = fP4;
+  Float_t oldfC00 = fC00;
+  Float_t oldfC10 = fC10;
+  Float_t oldfC11 = fC11;
+  Float_t oldfC20 = fC20;
+  Float_t oldfC21 = fC21;
+  Float_t oldfC22 = fC22;
+  Float_t oldfC30 = fC30;
+  Float_t oldfC31 = fC31;
+  Float_t oldfC32 = fC32;
+  Float_t oldfC33 = fC33;
+  Float_t oldfC40 = fC40;
+  Float_t oldfC41 = fC41;
+  Float_t oldfC42 = fC42;
+  Float_t oldfC43 = fC43;
+  Float_t oldfC44 = fC44;
+
+  if (TMath::Abs(fP4*Xnew - fP2) >= 0.9) // What's this??
+    {
+      //cout << "Propagation failed! Stiff track. " << pos << endl;
+      return 0;
+    }
+
+  if (TMath::Abs(dx) > 5) return 0;
+
+  /*if (TMath::Abs(fP4*Xold - fP2) >= 0.9) // What's this??
+    {
+      return 0;
+      }*/
+
+  // R must be approx. of the radius of the circle based on
+  // the old and new spacepoint. What then is Cod and Cnew??
+  // C has something to do with curvature at least.
+  
+  Float_t Cold = fP4*Xold - fP2;
+
+  //if (TMath::Abs(Cold) >= 1.0) return 0;
+  /*{
+    cout << "Cold " << endl << fP4*Xnew - fP2 << endl;
+    return 0;
+    }*/
+  Float_t Rold = TMath::Sqrt(1 - Cold*Cold);
+  Float_t Cnew = fP4*Xnew - fP2; 
+
+  //if (TMath::Abs(Cnew) >= 1.0) Cnew = 0.9;
+    /*{
+    cout << "Cnew " << endl;
+    return 0;
+    }*/
+  Float_t Rnew = TMath::Sqrt(1 - Cnew*Cnew);
+  //if (Rold < 0.9) 
+  //cout << "Cold = " << Cold << " , Rold = " << Rold << " , Cnew = " << Cnew << " , Rnew = " << Rnew << endl;  
+  // Prediction of the y- and z- coordinate in the next plane
+  fP0 += dx*(Cold+Cnew)/(Rold+Rnew);
+  fP1 += dx*(Cold+Cnew)/(Cold*Rnew + Cnew*Rold)*fP3;
+  //cout << "Old point " << Yold << " " << Zold << endl;
+  //cout << "Propagate " << fP0 << " " << fP1 << endl;
+  //cout << "Measured  " << points[pos].fY << " " << points[pos].fZ << endl;
+
+  fX = Xnew;
+
+  // f = F - 1 //What is this??
+  // Must be the f-matrix for the prediction, as in eq 1 in ALICE Kalman paper
+  Float_t RR = Rold + Rnew;
+  Float_t CC = Cold + Cnew;
+  Float_t XX = Xold + Xnew;
+
+  Float_t f02 = -dx*(2*RR + CC*(Cold/Rold + Cnew/Rnew))/(RR*RR);
+  Float_t f04 = dx*(RR*XX + CC*(Cold*Xold/Rold + Cnew*Xnew/Rnew))/(RR*RR);
+  Float_t CR = Cold*Rnew + Cnew*Rold;
+  Float_t f12 = -dx*fP3*(2*CR + CC*(Cnew*(Cold/Rold)-Rold + Cold*(Cnew/Rnew)-Rnew))/(CR*CR);
+  Float_t f13 = dx*CC/CR;
+  Float_t f14 = dx*fP3*(CR*XX-CC*(Rold*Xnew-Cnew*Cold*Xold/Rold + Rnew*Xold-Cold*Cnew*Xnew/Rnew))/(CR*CR);
+
+  // b = C*ft // This?
+  Float_t b00=f02*fC20 + f04*fC40;
+  Float_t b01=f12*fC20 + f14*fC40 + f13*fC30;
+  Float_t b10=f02*fC21 + f04*fC41;
+  Float_t b11=f12*fC21 + f14*fC41 + f13*fC31;
+  Float_t b20=f02*fC22 + f04*fC42;
+  Float_t b21=f12*fC22 + f14*fC42 + f13*fC32;
+  Float_t b30=f02*fC32 + f04*fC43;
+  Float_t b31=f12*fC32 + f14*fC43 + f13*fC33;
+  Float_t b40=f02*fC42 + f04*fC44;
+  Float_t b41=f12*fC42 + f14*fC44 + f13*fC43;
+
+  //a = f*b = f*C*ft
+  Float_t a00 = f02*b20 + f04*b40;
+  Float_t a01 = f02*b21 + f04*b41;
+  Float_t a11 = f12*b21 + f14*b41+f13*b31;
+
+  //F*C*Ft = C + (a + b + bt) /This is the covariance matrix, the samll t
+  // means transform. Then F must be df/dx
+  fC00 += a00 + 2*b00;
+  fC10 += a01 + b01 + b10;
+  fC20 += b20;
+  fC30 += b30;
+  fC40 += b40;
+  fC11 += a11 + 2*b11;
+  fC21 += b21;
+  fC31 += b31;
+  fC41 += b41;
+
+  // Update the track parameters with the measured values of the new point
+  Int_t value = UpdateOfflineTrack(x,y,z,ey,ez);
+  
+  if (value == 0)
+    {
+      fP0 = Yold;
+      fP1 = Zold;
+      fP2 = par2old;
+      fP3 = par3old;
+      fP4 = par4old;
+      fC00 = oldfC00;
+      fC10 = oldfC10;
+      fC11 = oldfC11;
+      fC20 = oldfC20;
+      fC21 = oldfC21;
+      fC22 = oldfC22;
+      fC30 = oldfC30;
+      fC31 = oldfC31;
+      fC32 = oldfC32;
+      fC33 = oldfC33;
+      fC40 = oldfC40;
+      fC41 = oldfC41;
+      fC42 = oldfC42;
+      fC43 = oldfC43;
+      fC44 = oldfC44;
+      
+      return value;
+      }
+  
+  else
+    return value;
+  //return 1;
+  //return UpdateTrack(points, pos);
+
+}
+
+Int_t AliL3KalmanTrack::UpdateOfflineTrack(Double_t x, Double_t y, Double_t z, Double_t ey, Double_t ez)
+{
+  // Update the track parameters with the measured values
+  //fX = x;
+
+  // The errors from the measurement of the spacepoint
+  Float_t sigmaY2 = ey;
+  Float_t sigmaZ2 = ez;
+  Float_t sigmaYZ = 0;
+  sigmaY2 += fC00;
+  sigmaZ2 += fC11;
+  sigmaYZ += fC10;
+
+  Float_t det = sigmaY2*sigmaZ2 - sigmaYZ*sigmaYZ;
+  Float_t tmp = sigmaY2;
+  sigmaY2 = sigmaZ2/det;
+  sigmaZ2 = tmp/det;
+  sigmaYZ = -sigmaYZ/det;
+  // Kalman gain matrix
+  Float_t k00 = fC00*sigmaY2 + fC10*sigmaYZ, k01 = fC00*sigmaYZ + fC10*sigmaZ2;
+  Float_t k10 = fC10*sigmaY2 + fC11*sigmaYZ, k11 = fC10*sigmaYZ + fC11*sigmaZ2;
+  Float_t k20 = fC20*sigmaY2 + fC21*sigmaYZ, k21 = fC20*sigmaYZ + fC21*sigmaZ2;
+  Float_t k30 = fC30*sigmaY2 + fC31*sigmaYZ, k31 = fC30*sigmaYZ + fC31*sigmaZ2;
+  Float_t k40 = fC40*sigmaY2 + fC41*sigmaYZ, k41 = fC40*sigmaYZ + fC41*sigmaZ2;
+  //cout << "x = " << fX << endl;
+  // Deviation between the predicted and measured values of y and z
+  Float_t dy = y-fP0; //cout << "dy = " << dy;
+  Float_t dz = z-fP1; //cout << ", dz = " << dz << endl;
+  //cout << "Measured " << points[pos].fY << " " << points[pos].fZ << endl;
+  
+  // Prediction of fP2 and fP4
+  Float_t cur = fP4 + k40*dy + k41*dz;
+  Float_t eta = fP2 + k20*dy + k21*dz;
+  
+  if (TMath::Abs(cur*fX - eta) >= 0.9)
+    {
+      //cout << "Update failed! Stiff track. " << pos << endl;
+      return 0;
+    }
+
+  // Filtered state vector
+  fP0 += k00*dy + k01*dz;
+  fP1 += k10*dy + k11*dz;
+  fP2 = eta;
+  fP3 += k30*dy + k31*dz; //cout << "update " << fP3 << endl;
+  fP4 = cur;
+  //cout << "AliL3KalmanTrack::Update, " << fP0 << " " << fP1 << " " << fP2 << " " << fP3 << " " << fP4 << endl;
+  //cout << "Measured, " << points[pos].fY << " " << points[pos].fZ << endl;
+  
+  Float_t c01 = fC10;
+  Float_t c02 = fC20;
+  Float_t c03 = fC30;
+  Float_t c04 = fC40;
+  Float_t c12 = fC21;
+  Float_t c13 = fC31;
+  Float_t c14 = fC41;
+  
+  // Filtered covariance matrix
+  fC00 -= k00*fC00 + k01*fC10; fC10 -= k00*c01 + k01*fC11;
+  fC20 -= k00*c02 + k01*c12;   fC30 -= k00*c03 + k01*c13;
+  fC40 -= k00*c04 + k01*c14;
+
+  fC11 -= k10*c01 + k11*fC11;
+  fC21 -= k10*c02 + k11*c12;   fC31 -= k10*c03 + k11*c13;
+  fC41 -= k10*c04 + k11*c14;
+
+  fC22 -= k20*c02 + k21*c12;   fC32 -= k20*c03 + k21*c13;
+  fC42 -= k20*c04 + k21*c14;
+
+  fC33 -= k30*c03 + k31*c13;
+  fC43 -= k40*c03 + k41*c13;
+
+  fC44 -= k40*c04 + k41*c14;
+
+  //cout << "AliL3KalmanTrack::Update, error " << fC00 << " " << fC11 << " " << fC22 << " " << fC33 << " " << fC44 << endl;
+
+  sigmaY2 = sigmaY2*det;
+  sigmaZ2 = sigmaZ2*det;
+  sigmaYZ = sigmaYZ*det;
+  //cout << "AliL3KalmanTrack::Update, Chi2 = " << GetChisq() << endl;
+
+  fChisq = GetChisq() + GetChisqIncrementOfflineTrack(y,z,ey,ez);
+  //cout << "fChisq = " << fChisq << endl;
+
+  // Must at some point make an cut on chisq. Here?
+  //if (fChisq > fMaxChi2) return 0;
+  
+  return 1;
+  
+}
+
+Float_t AliL3KalmanTrack::GetChisqIncrementOfflineTrack(Double_t y, Double_t z, Double_t ey, Double_t ez)
+{
+  // This function calculates a predicted chi2 increment.
+  //-----------------------------------------------------------------
+  Float_t r00=ey, r01=0., r11=ez;
+  r00+=fC00; r01+=fC10; r11+=fC11;
+
+  Double_t det=r00*r11 - r01*r01;
+  /*if (TMath::Abs(det) < 1.e-10) {
+    Int_t n=GetNumberOfClusters();
+    if (n>4) cerr<<n<<" AliKalmanTrack warning: Singular matrix !\n";
+    return 1e10;
+    }*/
+  Double_t tmp=r00; r00=r11; r11=tmp; r01=-r01;
+  
+  Double_t dy=y - fP0, dz=z - fP1;
+  //cout << "AliTPCtrack::GetPredictedChi2, r00 = " << r00 << " r11 = " << r11 << " ro1 = " << r01 << " dy = " << dy << " dz = " << dz << endl;
+  return (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz)/det;
+}
index 70d64d7..fbd722c 100644 (file)
 
 #include "AliL3RootTypes.h"
 #include "AliL3Track.h"
+
+// includes for offline comparison, will be removed
+#include "AliTPCtrack.h"
+// includes for offline comparison, will be removed
+
 class AliL3SpacePointData;
 
 class AliL3KalmanTrack : public AliL3Track {
@@ -54,10 +59,9 @@ class AliL3KalmanTrack : public AliL3Track {
 
   AliL3KalmanTrack();
   virtual ~AliL3KalmanTrack();
-  Int_t Init(AliL3Track *track, AliL3SpacePointData *points, UInt_t pos);
-  Int_t MakeTrackSeed(AliL3SpacePointData *points1, UInt_t pos1, AliL3SpacePointData *points2, UInt_t pos2, AliL3SpacePointData *points3, UInt_t pos3);
-  Int_t Propagate(AliL3SpacePointData *points, UInt_t pos);
-  Int_t UpdateTrack(AliL3SpacePointData *points, UInt_t pos);
+  Int_t Init(AliL3Track *track, AliL3SpacePointData *points, UInt_t pos,Int_t slice);
+  Int_t Propagate(AliL3SpacePointData *points, UInt_t pos, Int_t slice);
+  Int_t UpdateTrack(AliL3SpacePointData *points, UInt_t pos, Int_t slice);
   Int_t UpdateTrackII(AliL3SpacePointData *points, UInt_t pos);
   void AddTrack();
   //  Float_t GetStateVector(Float_t xx[5]) const {
@@ -91,6 +95,7 @@ class AliL3KalmanTrack : public AliL3Track {
   Float_t GetC13() {return fC43;}
   Float_t GetC14() {return fC44;}
 
+  void SetX(Float_t f) {fX = f;}
   void SetX0(Float_t f) {fP0 = f;}
   void SetX1(Float_t f) {fP1 = f;}
   void SetX2(Float_t f) {fP2 = f;}
@@ -132,7 +137,8 @@ class AliL3KalmanTrack : public AliL3Track {
   fC44 = f[14];}
   void SetChisq(Float_t f) {fChisq = f;}
   void SetMaxChi2(Float_t f) {fMaxChi2 = f;}
-  Float_t GetChisqIncrement(AliL3SpacePointData *points, UInt_t pos);
+  //Float_t GetChisqIncrement(AliL3SpacePointData *points, UInt_t pos);
+  Float_t GetChisqIncrement(Float_t y, Float_t error_y, Float_t z, Float_t error_z);
   Float_t f2(Float_t x1,Float_t y1, Float_t x2,Float_t y2, Float_t x3,Float_t y3);
   Float_t f3(Float_t x1,Float_t y1, Float_t x2,Float_t y2, Float_t z1,Float_t z2);
   Float_t f4(Float_t x1,Float_t y1, Float_t x2,Float_t y2, Float_t x3,Float_t y3);
@@ -140,6 +146,10 @@ class AliL3KalmanTrack : public AliL3Track {
 
   Int_t GetNHits() const {return fNHits;}
   void SetNHits(Int_t f) {fNHits = f;}
+
+  Int_t PropagateOfflineTrack(Double_t x, Double_t y, Double_t z, Double_t ey, Double_t ez);
+  Int_t UpdateOfflineTrack(Double_t x, Double_t y, Double_t z, Double_t ey, Double_t ez);
+  Float_t GetChisqIncrementOfflineTrack(Double_t y, Double_t z, Double_t ey, Double_t ez);
 };
 
 #endif
index c29d3a5..6dd0e4f 100644 (file)
@@ -19,7 +19,7 @@
 #endif
 #include "AliL3DDLDataFileHandler.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 50f4167..ae0a696 100644 (file)
@@ -9,6 +9,10 @@
 
 #include "AliL3DDLRawReaderFile.h"
 
+#if __GNUC__ == 3
+using namespace std;
+#endif
+
 /**************************************************************************
  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  *                                                                        *
 </pre>
 */
 
-
 ClassImp(AliL3DDLRawReaderFile)
 
-
 AliL3DDLRawReaderFile::AliL3DDLRawReaderFile(const Char_t* name, Bool_t addnum)
 {
   // create an object to read digits from the given input file(s)
index 3f46b6d..da776a2 100644 (file)
@@ -13,7 +13,7 @@
 #include "AliL3Transform.h"
 #include "AliL3DataHandler.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 3fbe309..316c154 100644 (file)
@@ -9,7 +9,7 @@
 
 #include "AliL3StandardIncludes.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index b902e44..a306b2c 100644 (file)
@@ -8,7 +8,7 @@
 #include "AliL3Transform.h"
 #include "AliL3TPCMapping.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 80fc352..9a3df5a 100644 (file)
@@ -7,7 +7,7 @@
 
 #include "AliL3TransBit.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 9da3d39..5390dab 100644 (file)
@@ -11,7 +11,7 @@
 //#define VHDLDEBUG
 #include "AliL3VHDLClusterFinder.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index d542844..4334179 100644 (file)
@@ -75,13 +75,9 @@ DEFSTR += -D$(ARCH) $(EXTRADEF)
 
 CXX          = g++
 LD           = $(CXX)
-GCCVERSION   = $(shell $(CXX) --version | head -n 1 | cut -d" " -f 3 | cut -d. -f 1 | cut -d" " -f1)
-
-CXXGCC3FLAGS = -DGCCVERSION=$(GCCVERSION)
 #CXXGCC3FLAGS += -pedantic
 #CXXGCC3FLAGS += -Wno-deprecated
 #CXXGCC3FLAGS += -Woverloaded-virtual
-
 CXXFLAGS       = -O2 -Wall -ggdb $(CXXGCC3FLAGS) $(EXTRACXXFLAGS) 
 LDFLAGS        = -O2 $(EXTRALDFLAGS) $(LIBS)
 STATICLDFLAGS  = -O2 $(PROFILEFLAGS) $(EXTRALDFLAGS) $(STATICLIBS)
index 00c0d85..14e530c 100644 (file)
@@ -15,7 +15,7 @@ BOMB THE COMPILE: USEPACKAGE=ALIROOT
 #include "AliL3DigitData.h"
 #include "AliL3Transform.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #else
 #include <stream.h>
index 0108ea3..f6b4aef 100644 (file)
@@ -10,7 +10,7 @@
 #include "AliL3StandardIncludes.h"
 #include "AliL3RootTypes.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 6eee667..cdf9235 100644 (file)
@@ -17,7 +17,7 @@
 #include <AliL3MemHandler.h>
 #include <AliL3HoughTransformerVhdl.h>
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
@@ -48,7 +48,7 @@ int main(Int_t argc,Char_t **argv)
   }
 
   AliL3Transform::Init(path);
-  //cerr << "Transform version: " << AliL3Transform::GetVersion() << endl;
+  cerr << "Transform version: " << AliL3Transform::GetVersion() << endl;
 
   AliL3HoughTransformerVhdl vtest(slice,patch,100,10);
   vtest.CreateHistograms(64,0.1,64,-30,30);
index 070cfaf..97114d6 100644 (file)
@@ -14,7 +14,7 @@
 #include "AliL3DigitData.h"
 #include "AliL3Transform.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #else
 #include <stream.h>
index bd85cb9..f793d36 100644 (file)
@@ -31,7 +31,7 @@
 #include <TGraphErrors.h>
 #endif
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 8129ab4..4984e02 100644 (file)
@@ -21,7 +21,7 @@
 #include "AliL3MemHandler.h"
 #include "AliL3SpacePointData.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #else
 #include <stream.h>
index c8bc13a..5e549ed 100644 (file)
@@ -10,7 +10,7 @@
 
 //Standalone program to run the track follower for benchmark tests.
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #else
 #include <stream.h>
index d31d6a6..58725f3 100644 (file)
@@ -15,7 +15,7 @@
 #include "AliL3AltroMemHandler.h"
 
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #else
 #include <stream.h>
index a28297a..8193e6e 100644 (file)
@@ -13,7 +13,7 @@
 #include "AliL3Logging.h"
 #include "AliL3Logger.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #else
 #include <stream.h>
index 5e4a765..be01eab 100644 (file)
@@ -12,7 +12,7 @@
 #include "AliL3SpacePointData.h"
 #include "AliL3MemHandler.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 2af01d6..05e12a7 100644 (file)
@@ -342,37 +342,13 @@ Int_t AliL3ConfMapFit::FitCircle()
 
   fTrack->SetCharge(q);
   
-//
-//    Get other track parameters
-//
-  Double_t x0, y0,phi0,r0,psi,pt ;
-  if ( fTrack->ComesFromMainVertex() == true ) 
-    {
-      //flag = 1 ; // primary track flag
-      x0   = fVertex->GetX() ;
-      y0   = fVertex->GetY() ;
-      phi0 = fVertex->GetPhi() ;
-      r0   = fVertex->GetR() ;
-      fTrack->SetPhi0(phi0);
-      fTrack->SetR0(r0);
-    } 
-  else 
-    {
-      //AliL3ConfMapPoint *lHit = (AliL3ConfMapPoint*)hits->Last();
-      AliL3ConfMapPoint *lHit = (AliL3ConfMapPoint*)fTrack->lastHit;
-      //flag =  0 ; // primary track flag
-      x0   =  lHit->GetX()  ;
-      y0   =  lHit->GetY()  ;
-      phi0 =  atan2(lHit->GetY(),lHit->GetX());
-      if ( phi0 < 0 ) phi0 += AliL3Transform::TwoPi();
-      r0   =  sqrt ( lHit->GetX() * lHit->GetX() + lHit->GetY() * lHit->GetY() )  ;
-      fTrack->SetPhi0(phi0);
-      fTrack->SetR0(r0);
-
-    }
   
   //Set the first point on the track to the space point coordinates of the innermost track
   //This will be updated to lie on the fit later on (AliL3Track::UpdateToFirstPoint).
+  Double_t x0,y0,psi,pt ;
+  AliL3ConfMapPoint *lHit = (AliL3ConfMapPoint*)fTrack->lastHit;
+  x0 = lHit->GetX();
+  y0 = lHit->GetY();
   fTrack->SetFirstPoint(x0,y0,0); //Z-value is set in FitLine
 
   psi  = (Double_t)atan2(bcent-y0,acent-x0) ;
@@ -418,7 +394,7 @@ Int_t AliL3ConfMapFit::FitLine ( )
   //TObjArray *hits = fTrack->GetHits();
   //Int_t num_of_hits = fTrack->GetNumberOfPoints();
 
-  if ( fTrack->ComesFromMainVertex() == true ) 
+  if (0)// fTrack->ComesFromMainVertex() == true ) 
     {
       dx = ((AliL3ConfMapPoint*)fTrack->firstHit)->GetX() - fVertex->GetX();
       dy = ((AliL3ConfMapPoint*)fTrack->firstHit)->GetY() - fVertex->GetY() ;
index 195748e..9d7973d 100644 (file)
@@ -182,10 +182,7 @@ void AliL3ConfMapTrack::Fill(AliL3Vertex *vertex,Double_t max_Dca)
       AliL3ConfMapPoint *fHit = (AliL3ConfMapPoint*)firstHit;
       SetLastPoint(fHit->GetX(),fHit->GetY(),fHit->GetZ());
       
-      if(AliLevel3::IsTracksAtFirstPoint())//Set the track fit parameter at the innermost most point
-       {
-         UpdateToFirstPoint();
-       }
+      UpdateToFirstPoint();
       
       delete fit;
     }
index 71cafb5..2e4b6e2 100644 (file)
@@ -30,7 +30,7 @@
 #include "AliL3SpacePointData.h"
 #include "AliL3MemHandler.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index a81e0a4..76a53bf 100644 (file)
@@ -20,7 +20,7 @@
 #include <AliComplexCluster.h>
 #include <AliStack.h>
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 #include <fstream>
 #include <iosfwd>
 #else
@@ -35,7 +35,7 @@
 #include "AliL3TrackArray.h"
 #include "AliL3Evaluate.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 04bedb3..f0a60a4 100644 (file)
@@ -27,7 +27,7 @@
 #include "AliL3TrackArray.h"
 #include "AliL3FileHandler.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index ee414ee..8aa6da8 100644 (file)
@@ -451,43 +451,20 @@ Int_t AliL3Fitter::FitCircle()
 
   fTrack->SetCharge(q);
   
-//
-//    Get other track parameters
-//
-  Double_t x0, y0,phi0,r0,psi,pt ;
-  if ( fVertexConstraint == kTRUE)
-    {
-      //flag = 1 ; // primary track flag
-      x0   = fVertex->GetX() ;
-      y0   = fVertex->GetY() ;
-      phi0 = fVertex->GetPhi() ;
-      r0   = fVertex->GetR() ;
-      fTrack->SetPhi0(phi0);
-      fTrack->SetR0(r0);
-    } 
-  else 
-    {
-      Int_t lastid=fTrack->GetNHits()-1;
-      UInt_t id = hitnum[lastid];
-      Int_t slice = (id>>25) & 0x7f;
-      Int_t patch = (id>>22) & 0x7;
-      UInt_t pos = id&0x3fffff;
-      AliL3SpacePointData *points = fClusters[slice][patch];
-      
-      //flag =  0 ; // primary track flag
-      x0   =  points[pos].fX;
-      y0   =  points[pos].fY;
-      phi0 =  atan2(points[pos].fY,points[pos].fX);
-      if ( phi0 < 0 ) phi0 += 2*AliL3Transform::Pi();
-      r0   =  sqrt ( points[pos].fX * points[pos].fX + points[pos].fY*points[pos].fY);
-      fTrack->SetPhi0(phi0);
-      fTrack->SetR0(r0);
-    }
-  
   //Set the first point on the track to the space point coordinates of the innermost track
   //This will be updated to lie on the fit later on (AliL3Track::UpdateToFirstPoint).
+  Double_t x0,y0,psi,pt ;
+  Int_t lastid=fTrack->GetNHits()-1;
+  UInt_t id = hitnum[lastid];
+  Int_t slice = (id>>25) & 0x7f;
+  Int_t patch = (id>>22) & 0x7;
+  UInt_t pos = id&0x3fffff;
+  AliL3SpacePointData *points = fClusters[slice][patch];
+  x0   =  points[pos].fX;
+  y0   =  points[pos].fY;
   fTrack->SetFirstPoint(x0,y0,0); //Z-value is set in FitLine
-
+  
+  //Set the remaining fit parameters
   psi  = (Double_t)atan2(bcent-y0,acent-x0) ;
   psi  = psi + q * 0.5F * AliL3Transform::Pi() ;
   if ( psi < 0 ) psi = psi + 2*AliL3Transform::Pi();
@@ -532,7 +509,7 @@ Int_t AliL3Fitter::FitLine ( )
   Double_t fS[(fTrack->GetNHits())];
   Double_t *fZWeight = new Double_t[fTrack->GetNHits()];
   UInt_t *hitnum = fTrack->GetHitNumbers();
-  if (fVertexConstraint==kTRUE)
+  if (0)//fVertexConstraint==kTRUE)
     {
       UInt_t id = hitnum[0];
       Int_t slice = (id>>25) & 0x7f;
index 8040e74..17bdb52 100644 (file)
@@ -33,7 +33,7 @@ class AliL3Logger{
   MLUCLogServer *so; //!
   MLUCLogServer *se; //!
   MLUCLogServer *sm; //!
-#if GCCVERSION == 3
+#if __GNUC__ == 3
   std::ofstream *of; //!
 #else  
   ofstream *of; //!
index 2b34a6c..2a71750 100644 (file)
@@ -25,7 +25,7 @@ class AliL3Log{
   enum TLogCmd { kEnd, kPrec, kHex, kDec };
 };
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 #define LOG( lvl, origin, keyword ) std::cerr<<"["<<origin<<": "<<keyword<<"] "
 #define ENDLOG std::endl
 #else
index d477b42..1eb3cd7 100644 (file)
@@ -13,7 +13,7 @@
 
 #include "AliL3RawDataFileHandler.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index e53176d..81b63bc 100644 (file)
@@ -16,7 +16,7 @@
 
 #else
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 #include <cstdio>
 #include <cmath>
 #else
index cdc4dcb..38fd91c 100644 (file)
@@ -3,7 +3,7 @@
 #ifndef ALIL3STANDARDINCLUDESH
 #define ALIL3STANDARDINCLUDESH
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 #include <iostream>
 #include <fstream>
 
@@ -35,7 +35,7 @@ eg. in inline functions defined in header files */
 #define STDCERR cerr
 #define STDENDL endl
 
-#endif //GCCVERSION
+#endif //__GNUC__
 
 #endif
 
index de440c5..90c895e 100644 (file)
@@ -11,7 +11,7 @@
 #include "AliL3Transform.h"
 #include "AliL3Vertex.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
index 5646fe1..5ddcd85 100644 (file)
@@ -199,6 +199,12 @@ void AliL3TrackArray::FillTracks(Int_t ntracks, AliL3TrackSegmentData* tr){
     track->SetFirstPoint(trs->fX,trs->fY,trs->fZ);
     track->SetLastPoint(trs->fLastX,trs->fLastY,trs->fLastZ);
     track->SetHits( trs->fNPoints, trs->fPointIDs );
+#ifdef ROWHOUGH
+    if(GetTrackType()=='h') {
+      ((AliL3HoughTrack *)track)->SetWeight(trs->fWeight);
+      track->SetMCid(trs->fTrackID);
+    }
+#endif
     UChar_t *tmpP = (UChar_t*)trs;
     tmpP += sizeof(AliL3TrackSegmentData)+trs->fNPoints*sizeof(UInt_t);
     trs = (AliL3TrackSegmentData*)tmpP;
@@ -230,6 +236,12 @@ void AliL3TrackArray::FillTracks(Int_t ntracks, AliL3TrackSegmentData* tr,Int_t
     UChar_t *tmpP = (UChar_t*)trs;
     tmpP += sizeof(AliL3TrackSegmentData)+trs->fNPoints*sizeof(UInt_t);
     trs = (AliL3TrackSegmentData*)tmpP;
+#ifdef ROWHOUGH
+    if(GetTrackType()=='h') {
+      ((AliL3HoughTrack *)track)->SetWeight(trs->fWeight);
+      track->SetMCid(trs->fTrackID);
+    }
+#endif
   }
 }
 
@@ -270,7 +282,12 @@ UInt_t AliL3TrackArray::WriteTracks(AliL3TrackSegmentData* tr){
     tP->fTgl = track->GetTgl();
     tP->fCharge = track->GetCharge();
     tP->fNPoints = track->GetNHits();
-
+#ifdef ROWHOUGH
+    if(GetTrackType()=='h') {
+      tP->fWeight = ((AliL3HoughTrack *)track)->GetWeight();
+      tP->fTrackID = track->GetMCid();
+    }
+#endif
     pP = (UInt_t*)track->GetHitNumbers();
     for (UInt_t j=0;j<tP->fNPoints;j++){
       tP->fPointIDs[j] = pP[j];
index 0adde8b..e18c489 100644 (file)
@@ -17,6 +17,10 @@ struct AliL3TrackSegmentData
        Double_t fPsi;
         Double_t fTgl;
         Int_t fCharge;
+#ifdef ROWHOUGH
+        UInt_t  fWeight;
+        Int_t  fTrackID;
+#endif
        UInt_t  fNPoints;
        UInt_t  fPointIDs[0];
     };
index 291ad66..9d67107 100644 (file)
@@ -24,7 +24,7 @@
 #include "AliL3Logging.h"
 #include "AliL3Transform.h"
 
-#if GCCVERSION == 3
+#if __GNUC__ == 3
 using namespace std;
 #endif
 
@@ -74,10 +74,8 @@ using namespace std;
 </pre>
 */
 
-
 ClassImp(AliL3Transform)
 
-
 const Double_t AliL3Transform::fAnodeWireSpacing = 0.25; //Taken from the TDR
 const Double_t AliL3Transform::fBFACT = 0.0029980;       //Conversion Factor
 const Double_t AliL3Transform::fPi  =   3.141592653589793;
@@ -1708,7 +1706,7 @@ void AliL3Transform::Global2HLT(Float_t *xyz,Int_t slice,Int_t slicerow)
 
 void AliL3Transform::PrintCompileOptions()
 {
-  cout << "Compiler (g++) version used: " << GCCVERSION << endl;
+  cout << "Compiler (g++) version used: " << __GNUC__ << endl;
 
 #ifdef no_root
   cout << "STANDALONE version: -Dno_root was given." << endl;
@@ -1739,9 +1737,22 @@ void AliL3Transform::PrintCompileOptions()
   cout << "NOT using Monte Carlo Info: -Ddo_mc was not given." << endl;
 #endif
 
+#ifdef ROWHOUGH
+  cout << "Using row transformer version: -DROWHOUGH was given." << endl;
+#else
+  cout << "NOT using row transformer version: -DROWHOUGH was not given." << endl;
+#endif
+
+#ifdef use_newio
+  cout << "Using NEWIO version: -Duse_newio was given." << endl;
+#else
+  cout << "NOT using NEWIO version: -Duse_newio was not given." << endl;
+#endif
+
 #ifdef use_logging
   cout << "Using logging classes (MLUC): -Duse_logging was given." << endl;
 #else
   cout << "NOT using logging classes (MLUC): -Duse_logging not was given." << endl;
 #endif
+
 }
index 3a54954..0407da1 100644 (file)
@@ -55,7 +55,6 @@
 
 ClassImp(AliLevel3)
 
-Bool_t AliLevel3::fSetTracks2FirstPoint = kTRUE;//Define track parameters at first point
 Bool_t AliLevel3::fDoVertexFit = kTRUE;//Include the vertex in the final track fit
 
 AliLevel3::AliLevel3()
@@ -567,8 +566,7 @@ void AliLevel3::FitGlobalTracks()
       AliL3Track *tr = tracks->GetCheckedTrack(i);
       if(!tr) continue;
       fitter->FitHelix(tr);
-      if(AliLevel3::IsTracksAtFirstPoint())
-       tr->UpdateToFirstPoint();
+      tr->UpdateToFirstPoint();
     }
   fBenchmark->Stop("Global track fitter");
   delete fitter;
index 0553bff..ad80b21 100644 (file)
@@ -68,9 +68,6 @@ class AliLevel3 : public TObject {
   Bool_t fUseBinary;
   Bool_t fWriteOut;
   
-  //Define whether track parameters should be given at first point on track (default)
-  //If not, the parameters will be given at the vertex.
-  static Bool_t fSetTracks2FirstPoint; 
   static Bool_t fDoVertexFit;
 
   Bool_t fClusterDeconv;
@@ -115,10 +112,7 @@ class AliLevel3 : public TObject {
   void DoRoi(Float_t e0=0.4,Float_t e1=0.5){fEta[0]=e0;fEta[1]=e1;fDoRoi=kTRUE;}
   void WriteFiles(Char_t *path="./"){fWriteOut = kTRUE; sprintf(fWriteOutPath,"%s",path);}
   
-  static void SetTracks2FirstPoint()   {fSetTracks2FirstPoint = kTRUE;}
-  static void SetTracks2Vertex()       {fSetTracks2FirstPoint = kFALSE;}
   static void SetVertexFit(Bool_t f)   {fDoVertexFit=f;}
-  static Bool_t IsTracksAtFirstPoint() {return fSetTracks2FirstPoint;}
   static Bool_t DoVertexFit()          {return fDoVertexFit;}
 
   ClassDef(AliLevel3,1) //Interface class for Level3-tracking
index 9b125fa..ff81e85 100644 (file)
@@ -3,6 +3,7 @@
 #include "AliITStrackV2.h"
 #include "AliL3Transform.h"
 #include <TVector3.h>
+#include <iostream.h>
 
 ClassImp(AliD0Trigger)
 
@@ -11,7 +12,7 @@ AliD0Trigger::AliD0Trigger()
   posTrack=0;
   negTrack=0;
 }
-AliD0Trigger::AliD0Trigger(double c[4],double b,double pv[3])
+AliD0Trigger::AliD0Trigger(double c[6],double b,double pv[3])
 {
   posTrack=0;
   negTrack=0;
@@ -19,6 +20,9 @@ AliD0Trigger::AliD0Trigger(double c[4],double b,double pv[3])
   cutV0high=c[1];
   cutInvMass=c[2]; 
   cutPointAngle=c[3];
+  cutd0d0=c[4];
+  cutCosThetaStar=c[5];
+  cutpTchild=c[6];
   Bfield=b;
   primaryVertex[0]=pv[0];
   primaryVertex[1]=pv[1];
@@ -72,74 +76,91 @@ bool AliD0Trigger::FindInvMass() {
 }
 bool AliD0Trigger::FindV0(){
 
-  bool goodV0=false;   
+  bool goodV0=false;
+
+  double Gxpos=0,Gypos=0,Gxneg=0,Gyneg=0;
+  double r1=0,r2=0,a1=0,a2=0,b1=0,b2=0;
+  double q=0,p=0,t=0,Fa=0,Fb=0,Fc=0;
+  double y1=0,y2=0,x1=0,x2=0,x3=0,x4=0;
+  double V01=0,V02=0,V03=0,V04=0;
+  double V0[4]={0,0,0,0};
+  bestV0[0]=0;  bestV0[1]=0;  bestV0[2]=0;
   
-  double Gxpos=posTrack->GetX()*cos(posTrack->GetAlpha())-posTrack->GetY()*sin(posTrack->GetAlpha());
-  double Gypos=posTrack->GetX()*sin(posTrack->GetAlpha())+posTrack->GetY()*cos(posTrack->GetAlpha());
-  double Gxneg=negTrack->GetX()*cos(negTrack->GetAlpha())-negTrack->GetY()*sin(negTrack->GetAlpha());
-  double Gyneg=negTrack->GetX()*sin(negTrack->GetAlpha())+negTrack->GetY()*cos(negTrack->GetAlpha());
+  Gxpos=posTrack->GetX()*cos(posTrack->GetAlpha())-posTrack->GetY()*sin(posTrack->GetAlpha());
+  Gypos=posTrack->GetX()*sin(posTrack->GetAlpha())+posTrack->GetY()*cos(posTrack->GetAlpha());
+  Gxneg=negTrack->GetX()*cos(negTrack->GetAlpha())-negTrack->GetY()*sin(negTrack->GetAlpha());
+  Gyneg=negTrack->GetX()*sin(negTrack->GetAlpha())+negTrack->GetY()*cos(negTrack->GetAlpha());
     
-  double r1=fabs(1/(AliL3Transform::GetBFact()*Bfield*posTrack->Get1Pt()));
-  double r2=fabs(1/(AliL3Transform::GetBFact()*Bfield*negTrack->Get1Pt()));
-  //double a1=posTrack->GetX()-(r1*cos(asin(posTrack->GetSnp())+TMath::PiOver2()));
-  //double a2=negTrack->GetX()-(r1*cos(asin(negTrack->GetSnp())-TMath::PiOver2()));
-  //double b1=posTrack->GetY()-(r1*sin(asin(posTrack->GetSnp())+TMath::PiOver2()));
-  //double b2=negTrack->GetY()-(r1*sin(asin(negTrack->GetSnp())-TMath::PiOver2()));
-  
-  double a1=Gxpos-(r1*cos(posTrack->GetAlpha()+asin(posTrack->GetSnp())+TMath::PiOver2()));
-  double a2=Gxneg-(r2*cos(negTrack->GetAlpha()+asin(negTrack->GetSnp())-TMath::PiOver2()));
-  double b1=Gypos-(r1*sin(posTrack->GetAlpha()+asin(posTrack->GetSnp())+TMath::PiOver2()));
-  double b2=Gyneg-(r2*sin(negTrack->GetAlpha()+asin(negTrack->GetSnp())-TMath::PiOver2()));
-
-  double a3=a1-a2;
-  double b3=b1-b2;
-  double r3=sqrt(pow(a3,2)+pow(b3,2));
-  
-  if(r3<r1+r2){
-    double q=a1-a2;
-    double p=2*(b1-b2);
-    double t=pow(r2,2)-pow(r2,2)+pow(b1,2)-pow(b2,2)-pow(q,2);
-    double Fa=pow(p,2)+(4*pow(q,2));
-    double Fb=-(2*t*p+8*pow(q,2)*b1);
-    double Fc=pow(t,2)-(4*pow(q,2)*pow(r1,2))+(4*pow(q,2)*pow(b1,2));
-    if(pow(Fb,2)-(4*Fa*Fc)>=0){
-      double y1=(-Fb+(sqrt(pow(Fb,2)-(4*Fa*Fc))))/(2*Fa);  //noe feil her. floating point
-      double y2=(-Fb-(sqrt(pow(Fb,2)-(4*Fa*Fc))))/(2*Fa);  //trolig negativ under rot 
-      
-      double x1=sqrt(pow(r1,2)-pow((y1-b1),2))+a1;    
-      double x2=sqrt(pow(r1,2)-pow((y2-b1),2))+a1;    
-      double x3=-sqrt(pow(r1,2)-pow((y1-b1),2))+a1;
-      double x4=-sqrt(pow(r1,2)-pow((y2-b1),2))+a1;
-      
-      double V01=sqrt(pow(x1,2)+pow(y1,2));
-      double V02=sqrt(pow(x2,2)+pow(y2,2));
-      double V03=sqrt(pow(x3,2)+pow(y1,2));
-      double V04=sqrt(pow(x4,2)+pow(y2,2));
-      
-      double V0[4]={V01,V02,V03,V04};
-      int index=0;
-      
-      double nearestV0=V0[0];
-      for(int i=1;i<4;i++){
-       if(nearestV0>V0[i]){
-         nearestV0=V0[i];
-         index=i;
-       }
+  r1=fabs(1/(AliL3Transform::GetBFact()*Bfield*posTrack->Get1Pt()));
+  r2=fabs(1/(AliL3Transform::GetBFact()*Bfield*negTrack->Get1Pt()));
+  
+  a1=Gxpos-(r1*cos(posTrack->GetAlpha()+asin(posTrack->GetSnp())+TMath::PiOver2()));
+  a2=Gxneg-(r2*cos(negTrack->GetAlpha()+asin(negTrack->GetSnp())-TMath::PiOver2()));
+  b1=Gypos-(r1*sin(posTrack->GetAlpha()+asin(posTrack->GetSnp())+TMath::PiOver2()));
+  b2=Gyneg-(r2*sin(negTrack->GetAlpha()+asin(negTrack->GetSnp())-TMath::PiOver2()));
+
+  //double a3=a1-a2;
+  //double b3=b1-b2;
+  //double r3=sqrt(pow(a3,2)+pow(b3,2));
+  
+  //if(r3<r1+r2){
+  q=a1-a2;
+  p=2*(b1-b2);
+  t=pow(r2,2)-pow(r2,2)+pow(b1,2)-pow(b2,2)-pow(q,2);
+  
+  Fa=pow(p,2)+(4*pow(q,2));
+  Fb=-(2*t*p+8*pow(q,2)*b1);
+  Fc=pow(t,2)-(4*pow(q,2)*pow(r1,2))+(4*pow(q,2)*pow(b1,2));
+  
+  if(pow(Fb,2)-(4*Fa*Fc)>=0){
+    y1=(-Fb+(sqrt(pow(Fb,2)-(4*Fa*Fc))))/(2*Fa);  //noe feil her. floating point
+    y2=(-Fb-(sqrt(pow(Fb,2)-(4*Fa*Fc))))/(2*Fa);  //trolig negativ under rot 
+    
+    x1=sqrt(pow(r1,2)-pow((y1-b1),2))+a1;    
+    x2=sqrt(pow(r1,2)-pow((y2-b1),2))+a1;    
+    x3=-sqrt(pow(r1,2)-pow((y1-b1),2))+a1;
+    x4=-sqrt(pow(r1,2)-pow((y2-b1),2))+a1;
+    
+    V01=sqrt(pow(x1-primaryVertex[0],2)+pow(y1-primaryVertex[1],2));
+    V02=sqrt(pow(x2-primaryVertex[0],2)+pow(y2-primaryVertex[1],2));
+    V03=sqrt(pow(x3-primaryVertex[0],2)+pow(y1-primaryVertex[1],2));  
+    V04=sqrt(pow(x4-primaryVertex[0],2)+pow(y2-primaryVertex[1],2));  
+    
+    V0[0]=V01;V0[1]=V02;V0[2]=V03;V0[3]=V04;
+
+    int index=0;
+    
+    double nearestV0=V0[0];
+
+    for(int i=0;i<4;i++){                        
+      if(nearestV0>V0[i]){
+       nearestV0=V0[i];
+       index=i;
       }
-      if(nearestV0<cutV0high && nearestV0>cutV0low){  
-       goodV0=true;
-       switch (index){
-       case 0 : bestV0[0]=x1;bestV0[1]=y1;bestV0[0]=0.;break;
-       case 1 : bestV0[0]=x2;bestV0[1]=y2;bestV0[0]=0.;break;
-       case 2 : bestV0[0]=x3;bestV0[1]=y1;bestV0[0]=0.;break;
-       case 3 : bestV0[0]=x4;bestV0[1]=y2;bestV0[0]=0.;break;
-       default : printf("Can't set V0");
+    }
+    cout<<"index: "<<index<<endl;
+    if(nearestV0<cutV0high){
+      if(nearestV0>cutV0low){  
+       goodV0=true;
+       if(index==0){
+         bestV0[0]=x1;   bestV0[1]=y1;   bestV0[2]=0.;
+       }
+       if(index==1){
+         bestV0[0]=x2;   bestV0[1]=y2;   bestV0[2]=0.;
+       }
+       if(index==2){
+         bestV0[0]=x3;   bestV0[1]=y1;   bestV0[2]=0.;
+       }
+       if(index==3){
+         bestV0[0]=x4;   bestV0[1]=y2;   bestV0[2]=0.;
        }
       }
     }
+    //}
   }
+  
+  cout<<"My V0: x: "<<bestV0[0]<<" y: "<<bestV0[1]<<" z: "<<bestV0[2]<<endl;
+
   return goodV0;       
 }
 void AliD0Trigger::FindMomentaAtVertex(){
@@ -205,3 +226,117 @@ Bool_t AliD0Trigger::PointingAngle(){
   else
     return false;
 }
+void AliD0Trigger::SetMomenta(double m[6])
+{
+  momenta[0] = m[0]; 
+  momenta[1] = m[1]; 
+  momenta[2] = m[2]; 
+  momenta[3] = m[3]; 
+  momenta[4] = m[4]; 
+  momenta[5] = m[5]; 
+}
+void AliD0Trigger::FindMomentaOffline()
+{
+  double ptP,alphaP,phiP,ptN,alphaN,phiN;
+  // momenta of the tracks at the vertex
+  ptP = 1./TMath::Abs(posTrack->Get1Pt());
+  alphaP = posTrack->GetAlpha();
+  phiP = alphaP+TMath::ASin(posTrack->GetSnp());
+  momenta[0] = ptP*TMath::Cos(phiP); 
+  momenta[1] = ptP*TMath::Sin(phiP);
+  momenta[2] = ptP*posTrack->GetTgl();
+  
+  ptN = 1./TMath::Abs(negTrack->Get1Pt());
+  alphaN = negTrack->GetAlpha();
+  phiN = alphaN+TMath::ASin(negTrack->GetSnp());
+  momenta[3] = ptN*TMath::Cos(phiN); 
+  momenta[4] = ptN*TMath::Sin(phiN);
+  momenta[5] = ptN*negTrack->GetTgl();  
+}
+bool AliD0Trigger::FindV0offline(double v[3])
+{
+
+  bool goodV0=false;
+  double r=0;
+  
+  bestV0[0]=v[0];
+  bestV0[1]=v[1];
+  bestV0[2]=v[2];
+
+  r=sqrt(pow(bestV0[0]-primaryVertex[0],2)+pow(bestV0[1]-primaryVertex[1],2)+pow(bestV0[2]-primaryVertex[2],2));
+
+  if(r<cutV0high){
+    if(r>cutV0low){
+      goodV0=true;
+    }
+  }
+
+  return goodV0;
+
+}
+bool AliD0Trigger::d0d0()
+{
+  bool goodd0=false;
+  double d00=0, d01=0;
+
+  d00 =  10000.*posTrack->GetD(primaryVertex[0],primaryVertex[1]);
+  d01 = -10000.*negTrack->GetD(primaryVertex[0],primaryVertex[1]);
+
+  if(d00*d01<cutd0d0){
+    goodd0=true;
+  }
+
+  return goodd0;
+}
+double AliD0Trigger::Energy()
+{
+  double kMD0 = 1.8645;  // D0  mass
+  return sqrt(P()*P()+kMD0*kMD0);
+}
+bool AliD0Trigger::CosThetaStar()
+{
+  bool goodtheta=false;
+  double kMD0 = 1.8645;  // D0  mass
+  double kMK  = 0.49368; // K+  mass
+  double kMPi = 0.13957; // Pi+ mass
+  double qL=0,qL2=0,ctsD0=0,ctsD0bar=0;
+
+  Double_t pStar = TMath::Sqrt(TMath::Power(kMD0*kMD0-kMK*kMK-kMPi*kMPi,2.)-4.*kMK*kMK*kMPi*kMPi)/(2.*kMD0);
+
+  Double_t beta = P()/Energy();
+  Double_t gamma = Energy()/kMD0;
+
+  TVector3 mom(momenta[0],momenta[1],momenta[2]);
+  TVector3 mom2(momenta[3],momenta[4],momenta[5]);
+  TVector3 momD(Px(),Py(),Pz());
+
+  qL = mom.Dot(momD)/momD.Mag();
+
+  ctsD0 = (qL/gamma-beta*TMath::Sqrt(pStar*pStar+kMK*kMK))/pStar;
+
+  qL2 = mom.Dot(momD)/momD.Mag();
+
+  ctsD0bar = (qL2/gamma-beta*TMath::Sqrt(pStar*pStar+kMK*kMK))/pStar;
+
+  if(TMath::Abs(ctsD0) > cutCosThetaStar){
+    if(TMath::Abs(ctsD0bar) > cutCosThetaStar){
+      goodtheta=true;
+    }
+  }
+  return goodtheta;
+}
+bool AliD0Trigger::pTchild()
+{
+  bool goodpT=false;
+  double pT1=0,pT2=0;
+  
+  pT1=sqrt(pow(momenta[0],2)+pow(momenta[1],2));
+  pT2=sqrt(pow(momenta[3],2)+pow(momenta[4],2));
+
+  if(pT1>cutpTchild){
+    if(pT2>cutpTchild){
+      goodpT=true;
+    }
+  }
+  return goodpT;
+}
index c9f3a90..576d9f9 100644 (file)
@@ -13,21 +13,33 @@ class AliD0Trigger {
   
   double momenta[6];
   double bestV0[3],primaryVertex[3];
-  double cutV0low, cutV0high, cutInvMass, cutPointAngle;
+  double cutV0low, cutV0high, cutInvMass, cutPointAngle, cutd0d0,cutCosThetaStar,cutpTchild;
   double Bfield;
 
  public:
   AliD0Trigger();
-  AliD0Trigger(double c[4],double Bfield,double pv[3]);
+  AliD0Trigger(double c[7],double Bfield,double pv[3]);
   AliD0Trigger(AliITStrackV2 * posT, AliITStrackV2 * negT);
   virtual ~AliD0Trigger();
 
   void SetTracks(AliITStrackV2 * posT, AliITStrackV2 * negT);
   bool FindInvMass();
   bool FindV0();
+  bool FindV0offline(double v[3]);
   void FindMomentaAtVertex();
+  void FindMomentaOffline();
   bool PointingAngle();
-  
+  void SetMomenta(double m[6]);
+  bool d0d0();
+  bool CosThetaStar();
+  double P(){return sqrt(Pt()*Pt()+Pz()*Pz());} 
+  double Pt(){return sqrt(Px()*Px()+Py()*Py());}
+  double Px(){return (momenta[0]+momenta[3]);}
+  double Py(){return (momenta[1]+momenta[4]);}
+  double Pz(){return (momenta[2]+momenta[5]);}
+  double Energy();
+  bool pTchild();
+
   ClassDef(AliD0Trigger,1) 
 
 };
index 54086a5..36f2ec0 100644 (file)
@@ -35,30 +35,17 @@ const Double_t kBz = 0.4;
 Double_t primaryvertex[3] = {0.,0.,0,};
 
 // sigle track cuts
-const Double_t kPtCut = 0.5;  // GeV/c
-const Double_t kd0Cut = 50.; // micron
+const Double_t kPtCut = 0.5;      // GeV/c
+const Double_t kd0Cut = 50.;      // micron
 const Double_t kd0CutHigh = 200.; // micron
 
 
 //cuts for combined tracks
-// cuts[0] = lowest V0 cut
-// cuts[1] = highest V0 cut
-// cuts[2] = inv. mass cut (diiferense)
-// cuts[3] = max value for pointing angle
-const Double_t cuts[4] = {0.005,
-                         0.02,
-                         0.5,
-                         0.95};
-
-// this function applies single track cuts
-Bool_t TrkCuts(const AliITStrackV2& trk);
-
-// this function creates TObjArrays with positive and negative tracks
-void   SelectTracks(TTree& itsTree,
-                      TObjArray& trksP,Int_t* itsEntryP,Int_t& nTrksP,
-                      TObjArray& trksN,Int_t* itsEntryN,Int_t& nTrksN);
-
-//void GetPrimaryVertex(int i,Char_t* path="./");
+const Double_t cuts[5] = {0.005,  // cuts[0] = lowest V0 cut  (cm)
+                         200,    // cuts[1] = highest V0 cut (cm) //0.02
+                         0.1,    // cuts[2] = inv. mass cut (diferense) (Gev/c)
+                         0.95,   // cuts[3] = max cosine value for pointing angle
+                          -5000}; // cuts[4] = max cosine value for pointing angle
 
 Int_t iTrkP,iTrkN,itsEntries;
 Char_t trksName[100];
@@ -66,7 +53,7 @@ Int_t nTrksP=0,nTrksN=0;
 Int_t nD0=0;
 int ev=0;
 
-void   RunD0Trigger(Int_t evFirst=0,Int_t evLast=1,Char_t* path="./") {
+void RunD0Trigger(Int_t evFirst=0,Int_t evLast=1,Char_t* path="./") {
 
   const Char_t *name="AliD0Trigger";
   cerr<<'\n'<<name<<"...\n";
diff --git a/HLT/trigger/RunD0offline.C b/HLT/trigger/RunD0offline.C
new file mode 100644 (file)
index 0000000..83d5f7d
--- /dev/null
@@ -0,0 +1,262 @@
+//#define __COMPILE__
+//#ifdef __COMPILE__
+#if !defined(__CINT__) || defined(__MAKECINT__)
+//-- --- standard headers------------- 
+#include <iostream.h>
+//--------Root headers ---------------
+#include <TSystem.h>
+#include <TFile.h>
+#include <TString.h>
+#include <TStopwatch.h>
+#include <TObject.h>
+#include <TVector3.h>
+#include <TTree.h>
+#include <TParticle.h>
+#include <TArray.h>
+//----- AliRoot headers ---------------
+#include "alles.h"
+#include "AliRun.h"
+#include "AliKalmanTrack.h"
+#include "AliITStrackV2.h"
+#include "AliHeader.h"
+#include "AliGenEventHeader.h"
+#include "AliV0vertex.h"
+#include "AliV0vertexer.h"
+#include "AliITSVertex.h"
+#include "AliITSVertexer.h"
+#include "AliITSVertexerTracks.h"
+#include "AliD0Trigger.h"
+#endif
+//-------------------------------------
+// field (T)
+const Double_t kBz = 0.4;
+
+// primary vertex
+Double_t primaryvertex[3] = {0.,0.,0,};
+
+//sec. vertex
+double v2[3]={0,0,0};
+
+// sigle track cuts
+const Double_t kPtCut = 0.5;  // GeV/c
+const Double_t kd0Cut = 50.; // micron
+const Double_t kd0CutHigh = 200.; // micron
+
+
+//cuts for combined tracks
+const Double_t cuts[7] = {0.005,     // cuts[0] = lowest V0 cut  (cm)
+                         0.015,     // cuts[1] = highest V0 cut (cm)
+                         0.05,      // cuts[2] = inv. mass cut (diferense) (Gev/c)
+                         0.95,      // cuts[3] = max cosine value for pointing angle
+                         -5000,     // cuts[4] = d0d0
+                         0.8,       // cuts[5] = costhetastar
+                         0.5};      // cuts[6] = ptchild
+//cut for distance of closest aprach
+double cutDCA=0.01;
+
+// this function applies single track cuts
+Bool_t TrkCuts(const AliITStrackV2& trk);
+
+// this function creates TObjArrays with positive and negative tracks
+void   SelectTracks(TTree& itsTree,
+                      TObjArray& trksP,Int_t* itsEntryP,Int_t& nTrksP,
+                      TObjArray& trksN,Int_t* itsEntryN,Int_t& nTrksN);
+
+//void GetPrimaryVertex(int i,Char_t* path="./");
+
+Int_t iTrkP,iTrkN,itsEntries;
+Char_t trksName[100];
+Int_t nTrksP=0,nTrksN=0;
+Int_t nD0=0;
+int ev=0;
+double mom[6];
+
+void RunD0offline(Int_t evFirst=0,Int_t evLast=1,Char_t* path="./") {
+
+  const Char_t *name="AliD0offline";
+  cerr<<'\n'<<name<<"...\n";
+  gBenchmark->Start(name); 
+
+  AliKalmanTrack::SetConvConst(100/0.299792458/kBz);
+
+  // Open file with ITS tracks
+  Char_t fnameTrack[1024];
+  sprintf(fnameTrack,"%s/AliITStracksV2.root",path);
+  TFile* itstrks = TFile::Open(fnameTrack);
+
+   // tracks from ITS
+  sprintf(trksName,"TreeT_ITS_%d",ev);
+  TTree *itsTree=(TTree*)itstrks->Get(trksName);
+  if(!itsTree) continue;
+  itsEntries = (Int_t)itsTree->GetEntries();
+  printf("+++\n+++ Number of tracks in ITS: %d\n+++\n\n",itsEntries);
+
+  //Getting primary Vertex
+  GetPrimaryVertex(0,path);
+
+  // call function which applies sigle track selection and
+  // separetes positives and negatives
+  TObjArray trksP(itsEntries/2);
+  Int_t *itsEntryP = new Int_t[itsEntries];
+  TObjArray trksN(itsEntries/2);
+  Int_t *itsEntryN = new Int_t[itsEntries];
+  SelectTracks(*itsTree,trksP,itsEntryP,nTrksP,trksN,itsEntryN,nTrksN); 
+
+  cout<<"#pos: "<<nTrksP<<endl;
+  cout<<"#neg: "<<nTrksN<<endl;
+
+  //the offline stuff
+  // define the cuts for vertexing
+  Double_t vtxcuts[]={33., // max. allowed chi2
+                     0.0, // min. allowed negative daughter's impact param 
+                     0.0, // min. allowed positive daughter's impact param 
+                     1.0, // max. allowed DCA between the daughter tracks
+                    -1.0, // min. allowed cosine of V0's pointing angle
+                     0.0, // min. radius of the fiducial volume
+                     2.9};// max. radius of the fiducial volume
+  
+  // create the AliV0vertexer object
+  AliV0vertexer *vertexer2 = new AliV0vertexer(vtxcuts);
+
+  AliD0Trigger * D0 = new AliD0Trigger(cuts,kBz,primaryvertex);
+
+  double ptP,alphaP,phiP,ptN,alphaN,phiN,dca;
+
+  for(iTrkP=0; iTrkP<nTrksP; iTrkP++) {
+    postrack = (AliITStrackV2*)trksP.At(iTrkP);
+    for(iTrkN=0; iTrkN<nTrksN; iTrkN++) {
+      negtrack = (AliITStrackV2*)trksN.At(iTrkN);
+      D0.SetTracks(postrack,negtrack);
+      //
+      // ----------- DCA MINIMIZATION ------------------
+      //
+      // find the DCA and propagate the tracks to the DCA 
+      double dca = vertexer2->PropagateToDCA(negtrack,postrack);
+  
+      if(dca<cutDCA){
+       // define the AliV0vertex object
+       AliV0vertex *vertex2 = new AliV0vertex(*negtrack,*postrack);
+       // get position of the vertex
+       vertex2->GetXYZ(v2[0],v2[1],v2[2]);
+       delete vertex2; 
+       if(D0.pTchild()){
+         if(D0.d0d0()){
+           if(D0.FindV0offline(v2)){
+             
+             // momenta of the tracks at the vertex
+             ptP = 1./TMath::Abs(postrack->Get1Pt());
+             alphaP = postrack->GetAlpha();
+             phiP = alphaP+TMath::ASin(postrack->GetSnp());
+             mom[0] = ptP*TMath::Cos(phiP); 
+             mom[1] = ptP*TMath::Sin(phiP);
+             mom[2] = ptP*postrack->GetTgl();
+             
+             ptN = 1./TMath::Abs(negtrack->Get1Pt());
+             alphaN = negtrack->GetAlpha();
+             phiN = alphaN+TMath::ASin(negtrack->GetSnp());
+             mom[3] = ptN*TMath::Cos(phiN); 
+             mom[4] = ptN*TMath::Sin(phiN);
+             mom[5] = ptN*negtrack->GetTgl();
+             
+             D0.SetMomenta(mom);
+             
+             if(D0.FindInvMass()){
+               if(D0.CosThetaStar()){
+                 if(D0.PointingAngle()){
+                   nD0++;
+                 }
+               }
+             }
+           } 
+         }
+       }
+      }
+    }
+  } 
+  cout<<"#D0: "<<nD0<<endl;
+  gBenchmark->Stop(name);
+  gBenchmark->Show(name);
+}
+//___________________________________________________________________________
+void   SelectTracks(TTree& itsTree,
+                    TObjArray& trksP,Int_t* itsEntryP,Int_t& nTrksP,
+                    TObjArray& trksN,Int_t* itsEntryN,Int_t& nTrksN) {
+//
+// this function creates two TObjArrays with positive and negative tracks
+//
+  nTrksP=0,nTrksN=0;
+
+  Int_t entr = (Int_t)itsTree.GetEntries();
+
+  // trasfer tracks from tree to arrays
+  for(Int_t i=0; i<entr; i++) {
+
+    AliITStrackV2 *itstrack = new AliITStrackV2; 
+    itsTree.SetBranchAddress("tracks",&itstrack);
+
+    itsTree.GetEvent(i);
+
+    // single track selection
+    if(!TrkCuts(*itstrack)) { delete itstrack; continue; }
+
+    if(itstrack->Get1Pt()>0.) { // negative track
+      trksN.AddLast(itstrack);
+      itsEntryN[nTrksN] = i;
+      nTrksN++;
+    } else {                    // positive track
+      trksP.AddLast(itstrack);
+      itsEntryP[nTrksP] = i;
+      nTrksP++;
+    }
+
+  }
+
+  return;
+}
+//____________________________________________________________________________
+Bool_t TrkCuts(const AliITStrackV2& trk) {
+// 
+// this function tells if track passes some kinematical cuts  
+//
+  if(TMath::Abs(1./trk.Get1Pt()) < kPtCut)                return kFALSE;
+  if(TMath::Abs(10000.*trk.GetD(primaryvertex[0],primaryvertex[1])) < kd0Cut) return kFALSE;
+  if(TMath::Abs(10000.*trk.GetD(primaryvertex[0],primaryvertex[1])) > kd0CutHigh) return kFALSE;
+
+  return kTRUE;
+}
+//____________________________________________________________________________
+void GetPrimaryVertex(int i,Char_t* path="./") {
+
+  int event=i;
+
+  Char_t falice[1024];
+  sprintf(falice,"%s/galice.root",path);
+  TFile * galice = new TFile(falice);
+  
+  TDirectory * curdir;  
+
+  Char_t vname[20];
+  galice->cd();
+  
+  sprintf(vname,"Vertex_%d",event);
+  TArrayF o = 0;
+  o.Set(3);
+  AliHeader * header = 0;
+  
+  TTree * treeE = (TTree*)gDirectory->Get("TE");
+  treeE->SetBranchAddress("Header",&header);
+  treeE->GetEntry(event);
+  AliGenEventHeader* genHeader = header->GenEventHeader();
+  if(genHeader){
+    genHeader->PrimaryVertex(o);
+    primaryvertex[0] = (Double_t)o[0];
+    primaryvertex[1] = (Double_t)o[1];
+    primaryvertex[2] = (Double_t)o[2];
+  }
+  else{
+    printf("Can't find Header");
+  }
+  delete header;
+  delete galice;
+}