]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/AliMUONClusterFinderVS.cxx
Cleaning of the code :
[u/mrichter/AliRoot.git] / MUON / AliMUONClusterFinderVS.cxx
index a60a8ee22a88a84fb88c93cb3eb0454c4048c31c..d5173d6d13a6509d77234c04ea49f4d7082ba549 100644 (file)
  **************************************************************************/
 /*
 $Log$
+Revision 1.8  2000/07/03 11:54:57  morsch
+AliMUONSegmentation and AliMUONHitMap have been replaced by AliSegmentation and AliHitMap in STEER
+The methods GetPadIxy and GetPadXxy of AliMUONSegmentation have changed name to GetPadI and GetPadC.
+
+Revision 1.7  2000/06/28 15:16:35  morsch
+(1) Client code adapted to new method signatures in AliMUONSegmentation (see comments there)
+to allow development of slat-muon chamber simulation and reconstruction code in the MUON
+framework. The changes should have no side effects (mostly dummy arguments).
+(2) Hit disintegration uses 3-dim hit coordinates to allow simulation
+of chambers with overlapping modules (MakePadHits, Disintegration).
+
+Revision 1.6  2000/06/28 12:19:18  morsch
+More consequent seperation of global input data services (AliMUONClusterInput singleton) and the
+cluster and hit reconstruction algorithms in AliMUONClusterFinderVS.
+AliMUONClusterFinderVS becomes the base class for clustering and hit reconstruction.
+It requires two cathode planes. Small modifications in the code will make it usable for
+one cathode plane and, hence, more general (for test beam data).
+AliMUONClusterFinder is now obsolete.
+
+Revision 1.5  2000/06/28 08:06:10  morsch
+Avoid global variables in AliMUONClusterFinderVS by seperating the input data for the fit from the
+algorithmic part of the class. Input data resides inside the AliMUONClusterInput singleton.
+It also naturally takes care of the TMinuit instance.
+
+Revision 1.4  2000/06/27 16:18:47  gosset
+Finally correct implementation of xm, ym, ixm, iym sizes
+when at least three local maxima on cathode 1 or on cathode 2
+
+Revision 1.3  2000/06/22 14:02:45  morsch
+Parameterised size of xm[], ym[], ixm[], iym[] correctly implemented (PH)
+Some HP scope problems corrected (PH)
+
 Revision 1.2  2000/06/15 07:58:48  morsch
 Code from MUON-dev joined
 
@@ -39,9 +71,8 @@ Revised and extended SplitByLocalMaxima method (Isabelle Chevrot):
 #include "AliMUONClusterFinderVS.h"
 #include "AliMUONDigit.h"
 #include "AliMUONRawCluster.h"
-#include "AliMUONSegmentation.h"
+#include "AliSegmentation.h"
 #include "AliMUONResponse.h"
-#include "AliMUONHitMap.h"
 #include "AliMUONHitMapA1.h"
 #include "AliRun.h"
 #include "AliMUON.h"
@@ -53,20 +84,12 @@ Revised and extended SplitByLocalMaxima method (Isabelle Chevrot):
 #include <TGraph.h> 
 #include <TPostScript.h> 
 #include <TMinuit.h> 
+#include <TF1.h>
+
 #include <stdio.h>
 #include <iostream.h>
 
 //_____________________________________________________________________
-static AliMUONSegmentation* fgSegmentation[2];
-static AliMUONResponse*     fgResponse;
-static Int_t                fgix[500][2];
-static Int_t                fgiy[500][2];
-static Float_t              fgCharge[500][2];
-static Int_t                fgNbins[2];
-static Int_t                fgFirst=kTRUE;
-static Int_t                fgChargeTot[2];
-static Float_t              fgQtot[2];
-static TMinuit*             fgMyMinuit ;
 // This function is minimized in the double-Mathieson fit
 void fcnS2(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag);
 void fcnS1(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag);
@@ -75,30 +98,12 @@ void fcnCombiS2(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t if
 
 ClassImp(AliMUONClusterFinderVS)
 
-    AliMUONClusterFinderVS::AliMUONClusterFinderVS
-(AliMUONSegmentation *seg1, AliMUONSegmentation *seg2,
- AliMUONResponse *response, 
- TClonesArray *digits1, TClonesArray *digits2, 
- Int_t chamber)   
-    :AliMUONClusterFinder(seg1, response, digits1, chamber)
-{
-// Constructor
-    fSegmentation2=seg2;
-    fDigits2=digits2;
-    fNdigits2 = fDigits2->GetEntriesFast();
-    fHitMap2=0;
-    fTrack[0]=fTrack[1]=-1;
-    
-}
-
     AliMUONClusterFinderVS::AliMUONClusterFinderVS()
-       :AliMUONClusterFinder()
 {
 // Default constructor
-    fSegmentation2=0;
-    fDigits2=0;
-    fNdigits2 = 0;
-    fHitMap2 = 0;
+    fInput=AliMUONClusterInput::Instance();
+    fHitMap[0] = 0;
+    fHitMap[1] = 0;
     fTrack[0]=fTrack[1]=-1;
 }
 
@@ -109,43 +114,6 @@ AliMUONClusterFinderVS::AliMUONClusterFinderVS(
     ;
 }
 
-void AliMUONClusterFinderVS::SetDigits(TClonesArray *MUONdigits1, TClonesArray *MUONdigits2) {
-// Set pointers to digit lists 
-
-    fDigits=MUONdigits1;
-    fNdigits = fDigits->GetEntriesFast();
-    fDigits2=MUONdigits2;
-    fNdigits2 = fDigits2->GetEntriesFast();
-}
-
-// Get Segmentation
-AliMUONSegmentation*  AliMUONClusterFinderVS::Segmentation(Int_t i)
-{
-// Return pointer to segmentation of cathode plane number 1 (i=0) or 2 (i=1)
-    return ((i==0)? fSegmentation : fSegmentation2);
-}
-
-// Get Number of Digits
-Int_t   AliMUONClusterFinderVS::NDigits(Int_t i)
-{
-// Return number of digits for cathode plane i+1
-    return ((i==0)? fNdigits : fNdigits2);
-}
-
-// Get Digits
-TClonesArray*  AliMUONClusterFinderVS::Digits(Int_t i)
-{
-// Return pointer to digits for cathode plane i+1
-    return ((i==0)? fDigits : fDigits2);
-}
-    
-
-AliMUONHitMap*   AliMUONClusterFinderVS::HitMap(Int_t i)
-{
-// Return pointer to  HitMap
-    return ((i==0)? fHitMap : fHitMap2);
-}
-
 void AliMUONClusterFinderVS::Decluster(AliMUONRawCluster *cluster)
 {
 // Decluster by local maxima
@@ -157,44 +125,33 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
 // Split complex cluster by local maxima 
 
     Int_t cath, i;
-    
+
+    fInput->SetCluster(c);
+
     fMul[0]=c->fMultiplicity[0];
     fMul[1]=c->fMultiplicity[1];
 
 //
 //  dump digit information into arrays
 //
-    fgSegmentation[0]=Segmentation(0);
-    fgSegmentation[1]=Segmentation(1);
-    fgResponse    =fResponse;
-    fgNbins[0]=fMul[0];
-    fgNbins[1]=fMul[1];
-    Float_t qtot;
+
+    Float_t qtot, zdum;
     
     for (cath=0; cath<2; cath++) {
        qtot=0;
        for (i=0; i<fMul[cath]; i++)
        {
            // pointer to digit
-           fDig[i][cath]=(AliMUONDigit*)
-               (Digits(cath)->UncheckedAt(c->fIndexMap[i][cath]));
+           fDig[i][cath]=fInput->Digit(cath, c->fIndexMap[i][cath]);
            // pad coordinates
            fIx[i][cath]= fDig[i][cath]->fPadX;
            fIy[i][cath]= fDig[i][cath]->fPadY;
            // pad charge
            fQ[i][cath] = fDig[i][cath]->fSignal;
            // pad centre coordinates
-           Segmentation(cath)->
-               GetPadCxy(fIx[i][cath], fIy[i][cath], fX[i][cath], fY[i][cath]);
-            // globals kUsed in fitting functions
-           fgix[i][cath]=fIx[i][cath];
-           fgiy[i][cath]=fIy[i][cath];
-           fgCharge[i][cath]=Float_t(fQ[i][cath]);
-           // total charge per cluster
-           qtot+=fgCharge[i][cath];
+           fInput->Segmentation(cath)->
+               GetPadC(fIx[i][cath], fIy[i][cath], fX[i][cath], fY[i][cath], zdum);
        } // loop over cluster digits
-       fgQtot[cath]=qtot;
-       fgChargeTot[cath]=Int_t(qtot);  
     }  // loop over cathodes
 
 
@@ -241,8 +198,8 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
        c->fY[1]=fYFit[0];
        c->fChi2[0]=chi2;
        c->fChi2[1]=chi2;
-       c->fX[0]=Segmentation(0)->GetAnod(c->fX[0]);
-       c->fX[1]=Segmentation(1)->GetAnod(c->fX[1]);
+       c->fX[0]=fInput->Segmentation(0)->GetAnod(c->fX[0]);
+       c->fX[1]=fInput->Segmentation(1)->GetAnod(c->fX[1]);
        
 // If reasonable chi^2 add result to the list of rawclusters
        //      if (chi2 < 50) {
@@ -334,12 +291,12 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
        for (ico=0; ico<4; ico++) {
            accepted[ico]=kFALSE;
 // cathode one: x-coordinate
-           isec=Segmentation(0)->Sector(ixm[ico][0], iym[ico][0]);
-           dpx=Segmentation(0)->Dpx(isec)/2.;
+           isec=fInput->Segmentation(0)->Sector(ixm[ico][0], iym[ico][0]);
+           dpx=fInput->Segmentation(0)->Dpx(isec)/2.;
            dx=TMath::Abs(xm[ico][0]-xm[ico][1]);
 // cathode two: y-coordinate
-           isec=Segmentation(1)->Sector(ixm[ico][1], iym[ico][1]);
-           dpy=Segmentation(1)->Dpy(isec)/2.;
+           isec=fInput->Segmentation(1)->Sector(ixm[ico][1], iym[ico][1]);
+           dpy=fInput->Segmentation(1)->Dpy(isec)/2.;
            dy=TMath::Abs(ym[ico][0]-ym[ico][1]);
 //         printf("\n %i %f %f %f %f \n", ico, ym[ico][0], ym[ico][1], dy, dpy );
            if ((dx <= dpx) && (dy <= dpy)) {
@@ -464,16 +421,16 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
                            AliMUONRawCluster cnew;
                            Int_t cath;    
                            for (cath=0; cath<2; cath++) {
-                             cnew.fX[cath]=Float_t(xm[ico][1]);
-                             cnew.fY[cath]=Float_t(ym[ico][0]);
-                             cnew.fMultiplicity[cath]=c->fMultiplicity[cath];
+                               cnew.fX[cath]=Float_t(xm[ico][1]);
+                               cnew.fY[cath]=Float_t(ym[ico][0]);
+                               cnew.fMultiplicity[cath]=c->fMultiplicity[cath];
                                for (i=0; i<fMul[cath]; i++) {
-                                 cnew.fIndexMap[i][cath]=c->fIndexMap[i][cath];
-                                 fgSegmentation[cath]->SetPad(fgix[i][cath], fgiy[i][cath]);
+                                   cnew.fIndexMap[i][cath]=c->fIndexMap[i][cath];
+                                   fInput->Segmentation(cath)->SetPad(fIx[i][cath], fIy[i][cath]);
                                }
-                             fprintf(stderr,"\nRawCluster %d cath %d\n",ico,cath);
-                             fprintf(stderr,"mult_av %d\n",c->fMultiplicity[cath]);
-                             FillCluster(&cnew,cath);
+                               fprintf(stderr,"\nRawCluster %d cath %d\n",ico,cath);
+                               fprintf(stderr,"mult_av %d\n",c->fMultiplicity[cath]);
+                               FillCluster(&cnew,cath);
                            } 
                            cnew.fClusterType=cnew.PhysicsContribution();
                            AddRawCluster(cnew);
@@ -483,90 +440,90 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
                }
            }
        }
-
+       
 // ******* iacc = 2 *******
 // Two combinations found between the 2 cathodes
        if (iacc==2) {
 
 // Was the same maximum taken twice
-               if ((accepted[0]&&accepted[1]) || (accepted[2]&&accepted[3])) {
-                   fprintf(stderr,"\n Maximum taken twice !!!\n");
+           if ((accepted[0]&&accepted[1]) || (accepted[2]&&accepted[3])) {
+               fprintf(stderr,"\n Maximum taken twice !!!\n");
 
 // Have a try !! with that 
-                   if (accepted[0]&&accepted[3]) {
-                       fXInit[0]=xm[0][1];
-                       fYInit[0]=ym[0][0];
-                       fXInit[1]=xm[1][1];
-                       fYInit[1]=ym[1][0];
-                   } else {
-                       fXInit[0]=xm[2][1];
-                       fYInit[0]=ym[2][0];
-                       fXInit[1]=xm[3][1];
-                       fYInit[1]=ym[3][0];
-                   }
-                   fprintf(stderr,"\n cas (2) CombiDoubleMathiesonFit(c)\n");
-                   chi2=CombiDoubleMathiesonFit(c);
+               if (accepted[0]&&accepted[3]) {
+                   fXInit[0]=xm[0][1];
+                   fYInit[0]=ym[0][0];
+                   fXInit[1]=xm[1][1];
+                   fYInit[1]=ym[1][0];
+               } else {
+                   fXInit[0]=xm[2][1];
+                   fYInit[0]=ym[2][0];
+                   fXInit[1]=xm[3][1];
+                   fYInit[1]=ym[3][0];
+               }
+               fprintf(stderr,"\n cas (2) CombiDoubleMathiesonFit(c)\n");
+               chi2=CombiDoubleMathiesonFit(c);
 //                 Int_t ndf = fgNbins[0]+fgNbins[1]-6;
 //                 Float_t prob = TMath::Prob(chi2,ndf);
 //                 prob2->Fill(prob);
 //                 chi2_2->Fill(chi2);
-                   Split(c);
-
-               } else {
+               Split(c);
+               
+           } else {
 // No ghosts ! No Problems ! -  Perform one fit only !
-                   if (accepted[0]&&accepted[3]) {
-                       fXInit[0]=xm[0][1];
-                       fYInit[0]=ym[0][0];
-                       fXInit[1]=xm[3][1];
-                       fYInit[1]=ym[3][0];
-                   } else {
-                       fXInit[0]=xm[1][1];
-                       fYInit[0]=ym[1][0];
-                       fXInit[1]=xm[2][1];
-                       fYInit[1]=ym[2][0];
-                   }
-                   fprintf(stderr,"\n cas (2) CombiDoubleMathiesonFit(c)\n");
-                   chi2=CombiDoubleMathiesonFit(c);
+               if (accepted[0]&&accepted[3]) {
+                   fXInit[0]=xm[0][1];
+                   fYInit[0]=ym[0][0];
+                   fXInit[1]=xm[3][1];
+                   fYInit[1]=ym[3][0];
+               } else {
+                   fXInit[0]=xm[1][1];
+                   fYInit[0]=ym[1][0];
+                   fXInit[1]=xm[2][1];
+                   fYInit[1]=ym[2][0];
+               }
+               fprintf(stderr,"\n cas (2) CombiDoubleMathiesonFit(c)\n");
+               chi2=CombiDoubleMathiesonFit(c);
 //                 Int_t ndf = fgNbins[0]+fgNbins[1]-6;
 //                 Float_t prob = TMath::Prob(chi2,ndf);
 //                 prob2->Fill(prob);
 //                 chi2_2->Fill(chi2);
-                   fprintf(stderr," chi2 %f\n",chi2);
-                   Split(c);
-               }
-
+               fprintf(stderr," chi2 %f\n",chi2);
+               Split(c);
+           }
+           
 // ******* iacc = 4 *******
 // Four combinations found between the 2 cathodes
 // Ghost !!
-           } else if (iacc==4) {
+       } else if (iacc==4) {
 // Perform fits for the two possibilities !!   
-               fXInit[0]=xm[0][1];
-               fYInit[0]=ym[0][0];
-               fXInit[1]=xm[3][1];
-               fYInit[1]=ym[3][0];
-               fprintf(stderr,"\n cas (2) CombiDoubleMathiesonFit(c)\n");
-               chi2=CombiDoubleMathiesonFit(c);
+           fXInit[0]=xm[0][1];
+           fYInit[0]=ym[0][0];
+           fXInit[1]=xm[3][1];
+           fYInit[1]=ym[3][0];
+           fprintf(stderr,"\n cas (2) CombiDoubleMathiesonFit(c)\n");
+           chi2=CombiDoubleMathiesonFit(c);
 //             Int_t ndf = fgNbins[0]+fgNbins[1]-6;
 //             Float_t prob = TMath::Prob(chi2,ndf);
 //             prob2->Fill(prob);
 //             chi2_2->Fill(chi2);
-               fprintf(stderr," chi2 %f\n",chi2);
-               Split(c);
-               fXInit[0]=xm[1][1];
-               fYInit[0]=ym[1][0];
-               fXInit[1]=xm[2][1];
-               fYInit[1]=ym[2][0];
-               fprintf(stderr,"\n cas (2) CombiDoubleMathiesonFit(c)\n");
-               chi2=CombiDoubleMathiesonFit(c);
+           fprintf(stderr," chi2 %f\n",chi2);
+           Split(c);
+           fXInit[0]=xm[1][1];
+           fYInit[0]=ym[1][0];
+           fXInit[1]=xm[2][1];
+           fYInit[1]=ym[2][0];
+           fprintf(stderr,"\n cas (2) CombiDoubleMathiesonFit(c)\n");
+           chi2=CombiDoubleMathiesonFit(c);
 //             ndf = fgNbins[0]+fgNbins[1]-6;
 //             prob = TMath::Prob(chi2,ndf);
 //             prob2->Fill(prob);
 //             chi2_2->Fill(chi2);
-               fprintf(stderr," chi2 %f\n",chi2);
-               Split(c);
-           }
+           fprintf(stderr," chi2 %f\n",chi2);
+           Split(c);
+       }
 
-       } else if (fNLocal[0]==2 &&  fNLocal[1]==1) {
+    } else if (fNLocal[0]==2 &&  fNLocal[1]==1) {
 //  ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 //  (3) Two local maxima on cathode 1 and one maximum on cathode 2 
 //  ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -580,16 +537,16 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
 //  0-0, 0-1, 1-0, 1-1 
         ico=0;
        for (im1=0; im1<2; im1++) {
-               xm[ico][0]=fX[fIndLocal[im1][0]][0];
-               ym[ico][0]=fY[fIndLocal[im1][0]][0];
-               xm[ico][1]=fX[fIndLocal[0][1]][1];
-               ym[ico][1]=fY[fIndLocal[0][1]][1];
-
-               ixm[ico][0]=fIx[fIndLocal[im1][0]][0];
-               iym[ico][0]=fIy[fIndLocal[im1][0]][0];
-               ixm[ico][1]=fIx[fIndLocal[0][1]][1];
-               iym[ico][1]=fIy[fIndLocal[0][1]][1];
-               ico++;
+           xm[ico][0]=fX[fIndLocal[im1][0]][0];
+           ym[ico][0]=fY[fIndLocal[im1][0]][0];
+           xm[ico][1]=fX[fIndLocal[0][1]][1];
+           ym[ico][1]=fY[fIndLocal[0][1]][1];
+           
+           ixm[ico][0]=fIx[fIndLocal[im1][0]][0];
+           iym[ico][0]=fIy[fIndLocal[im1][0]][0];
+           ixm[ico][1]=fIx[fIndLocal[0][1]][1];
+           iym[ico][1]=fIy[fIndLocal[0][1]][1];
+           ico++;
        }
 // ico = 0 : first local maximum on cathodes 1 and 2
 // ico = 1 : second local maximum on cathode 1 and first on cathode 2
@@ -602,11 +559,11 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
        
        for (ico=0; ico<2; ico++) {
            accepted[ico]=kFALSE;
-           isec=Segmentation(0)->Sector(ixm[ico][0], iym[ico][0]);
-           dpx=Segmentation(0)->Dpx(isec)/2.;
+           isec=fInput->Segmentation(0)->Sector(ixm[ico][0], iym[ico][0]);
+           dpx=fInput->Segmentation(0)->Dpx(isec)/2.;
            dx=TMath::Abs(xm[ico][0]-xm[ico][1]);
-           isec=Segmentation(1)->Sector(ixm[ico][1], iym[ico][1]);
-           dpy=Segmentation(1)->Dpy(isec)/2.;
+           isec=fInput->Segmentation(1)->Sector(ixm[ico][1], iym[ico][1]);
+           dpy=fInput->Segmentation(1)->Dpy(isec)/2.;
            dy=TMath::Abs(ym[ico][0]-ym[ico][1]);
 //         printf("\n %i %f %f %f %f \n", ico, ym[ico][0], ym[ico][1], dy, dpy );
            if ((dx <= dpx) && (dy <= dpy)) {
@@ -618,10 +575,10 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
                accepted[ico]=kFALSE;
            }
        }
-
+       
        Float_t chi21 = 100;
        Float_t chi22 = 100;
-
+       
        if (accepted[0]) {
            fXInit[0]=xm[0][1];
            fYInit[0]=ym[0][0];
@@ -659,8 +616,8 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
                        cnew.fY[cath]=Float_t(ym[ico][0]);
                        cnew.fMultiplicity[cath]=c->fMultiplicity[cath];
                        for (i=0; i<fMul[cath]; i++) {
-                         cnew.fIndexMap[i][cath]=c->fIndexMap[i][cath];
-                         fgSegmentation[cath]->SetPad(fgix[i][cath], fgiy[i][cath]);
+                           cnew.fIndexMap[i][cath]=c->fIndexMap[i][cath];
+                           fInput->Segmentation(cath)->SetPad(fIx[i][cath], fIy[i][cath]);
                        }
                        fprintf(stderr,"\nRawCluster %d cath %d\n",ico,cath);
                        fprintf(stderr,"mult_av %d\n",c->fMultiplicity[cath]);
@@ -672,12 +629,12 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
                }
            }
        }
-          
+       
 //  +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 //  (3') One local maximum on cathode 1 and two maxima on cathode 2 
 //  +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
     } else if (fNLocal[0]==1 && fNLocal[1]==2) {
-
+       
        Float_t xm[4][2], ym[4][2];
        Float_t dpx, dpy, dx, dy;
        Int_t ixm[4][2], iym[4][2];
@@ -687,16 +644,16 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
 //  0-0, 0-1, 1-0, 1-1 
         ico=0;
        for (im1=0; im1<2; im1++) {
-               xm[ico][0]=fX[fIndLocal[0][0]][0];
-               ym[ico][0]=fY[fIndLocal[0][0]][0];
-               xm[ico][1]=fX[fIndLocal[im1][1]][1];
-               ym[ico][1]=fY[fIndLocal[im1][1]][1];
-
-               ixm[ico][0]=fIx[fIndLocal[0][0]][0];
-               iym[ico][0]=fIy[fIndLocal[0][0]][0];
-               ixm[ico][1]=fIx[fIndLocal[im1][1]][1];
-               iym[ico][1]=fIy[fIndLocal[im1][1]][1];
-               ico++;
+           xm[ico][0]=fX[fIndLocal[0][0]][0];
+           ym[ico][0]=fY[fIndLocal[0][0]][0];
+           xm[ico][1]=fX[fIndLocal[im1][1]][1];
+           ym[ico][1]=fY[fIndLocal[im1][1]][1];
+           
+           ixm[ico][0]=fIx[fIndLocal[0][0]][0];
+           iym[ico][0]=fIy[fIndLocal[0][0]][0];
+           ixm[ico][1]=fIx[fIndLocal[im1][1]][1];
+           iym[ico][1]=fIy[fIndLocal[im1][1]][1];
+           ico++;
        }
 // ico = 0 : first local maximum on cathodes 1 and 2
 // ico = 1 : first local maximum on cathode 1 and second on cathode 2
@@ -709,11 +666,11 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
        
        for (ico=0; ico<2; ico++) {
            accepted[ico]=kFALSE;
-           isec=Segmentation(0)->Sector(ixm[ico][0], iym[ico][0]);
-           dpx=Segmentation(0)->Dpx(isec)/2.;
+           isec=fInput->Segmentation(0)->Sector(ixm[ico][0], iym[ico][0]);
+           dpx=fInput->Segmentation(0)->Dpx(isec)/2.;
            dx=TMath::Abs(xm[ico][0]-xm[ico][1]);
-           isec=Segmentation(1)->Sector(ixm[ico][1], iym[ico][1]);
-           dpy=Segmentation(1)->Dpy(isec)/2.;
+           isec=fInput->Segmentation(1)->Sector(ixm[ico][1], iym[ico][1]);
+           dpy=fInput->Segmentation(1)->Dpy(isec)/2.;
            dy=TMath::Abs(ym[ico][0]-ym[ico][1]);
 //         printf("\n %i %f %f %f %f \n", ico, ym[ico][0], ym[ico][1], dy, dpy );
            if ((dx <= dpx) && (dy <= dpy)) {
@@ -767,8 +724,8 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
                        cnew.fY[cath]=Float_t(ym[ico][0]);
                        cnew.fMultiplicity[cath]=c->fMultiplicity[cath];
                        for (i=0; i<fMul[cath]; i++) {
-                         cnew.fIndexMap[i][cath]=c->fIndexMap[i][cath];
-                         fgSegmentation[cath]->SetPad(fgix[i][cath], fgiy[i][cath]);
+                           cnew.fIndexMap[i][cath]=c->fIndexMap[i][cath];
+                           fInput->Segmentation(cath)->SetPad(fIx[i][cath], fIy[i][cath]);
                        }
                        fprintf(stderr,"\nRawCluster %d cath %d\n",ico,cath);
                        fprintf(stderr,"mult_av %d\n",c->fMultiplicity[cath]);
@@ -788,15 +745,15 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
        
        Int_t param = fNLocal[0]*fNLocal[1];
        Int_t ii;
-       
-       Float_t ** xm = new Float_t * [2];
-       for (ii=0; ii<2; ii++) xm[ii]=new Float_t [param];
-       Float_t ** ym = new Float_t * [2];
-       for (ii=0; ii<2; ii++) ym[ii]=new Float_t [param];
-       Int_t ** ixm = new Int_t * [2];
-       for (ii=0; ii<2; ii++) ixm[ii]=new Int_t [param];
-       Int_t ** iym = new Int_t * [2];
-       for (ii=0; ii<2; ii++) iym[ii]=new Int_t [param];
+
+       Float_t ** xm = new Float_t * [param];
+       for (ii=0; ii<param; ii++) xm[ii]=new Float_t [2];
+       Float_t ** ym = new Float_t * [param];
+       for (ii=0; ii<param; ii++) ym[ii]=new Float_t [2];
+       Int_t ** ixm = new Int_t * [param];
+       for (ii=0; ii<param; ii++) ixm[ii]=new Int_t [2];
+       Int_t ** iym = new Int_t * [param];
+       for (ii=0; ii<param; ii++) iym[ii]=new Int_t [2];
        
        Int_t isec, ico;
        Float_t dpx, dpy, dx, dy;
@@ -816,17 +773,17 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
                ico++;
            }
        }
-
+       
        Int_t nIco = ico;
-
+       
        fprintf(stderr,"nIco %d\n",nIco);
        for (ico=0; ico<nIco; ico++) {
            fprintf(stderr,"ico = %d\n",ico);
-           isec=Segmentation(0)->Sector(ixm[ico][0], iym[ico][0]);
-           dpx=Segmentation(0)->Dpx(isec)/2.;
+           isec=fInput->Segmentation(0)->Sector(ixm[ico][0], iym[ico][0]);
+           dpx=fInput->Segmentation(0)->Dpx(isec)/2.;
            dx=TMath::Abs(xm[ico][0]-xm[ico][1]);
-           isec=Segmentation(1)->Sector(ixm[ico][1], iym[ico][1]);
-           dpy=Segmentation(1)->Dpy(isec)/2.;
+           isec=fInput->Segmentation(1)->Sector(ixm[ico][1], iym[ico][1]);
+           dpy=fInput->Segmentation(1)->Dpy(isec)/2.;
            dy=TMath::Abs(ym[ico][0]-ym[ico][1]);
 
            fprintf(stderr,"dx %f dpx %f dy %f dpy %f\n",dx,dpx,dy,dpy);
@@ -840,8 +797,8 @@ void AliMUONClusterFinderVS::SplitByLocalMaxima(AliMUONRawCluster *c)
                    cnew.fY[cath]=Float_t(ym[ico][0]);
                    cnew.fMultiplicity[cath]=c->fMultiplicity[cath];
                    for (i=0; i<fMul[cath]; i++) {
-                     cnew.fIndexMap[i][cath]=c->fIndexMap[i][cath];
-                     fgSegmentation[cath]->SetPad(fgix[i][cath], fgiy[i][cath]);
+                       cnew.fIndexMap[i][cath]=c->fIndexMap[i][cath];
+                       fInput->Segmentation(cath)->SetPad(fIx[i][cath], fIy[i][cath]);
                    }
                    FillCluster(&cnew,cath);
                } 
@@ -878,23 +835,23 @@ void AliMUONClusterFinderVS::FindLocalMaxima(AliMUONRawCluster* c)
     }
 //  number of next neighbours and arrays to store them 
     Int_t nn;
-    Int_t x[kMaxNeighbours], y[kMaxNeighbours];
+    Int_t x[10], y[10];
 // loop over cathodes
     for (cath=0; cath<2; cath++) {
 // loop over cluster digits
        for (i=0; i<fMul[cath]; i++) {
 // get neighbours for that digit and assume that it is local maximum       
-           Segmentation(cath)->Neighbours(fIx[i][cath], fIy[i][cath], &nn, x, y);
+           fInput->Segmentation(cath)->Neighbours(fIx[i][cath], fIy[i][cath], &nn, x, y);
            isLocal[i][cath]=kTRUE;
-           Int_t isec= Segmentation(cath)->Sector(fIx[i][cath], fIy[i][cath]);
-           Float_t a0 = Segmentation(cath)->Dpx(isec)*Segmentation(cath)->Dpy(isec);
+           Int_t isec= fInput->Segmentation(cath)->Sector(fIx[i][cath], fIy[i][cath]);
+           Float_t a0 = fInput->Segmentation(cath)->Dpx(isec)*fInput->Segmentation(cath)->Dpy(isec);
 // loop over next neighbours, if at least one neighbour has higher charger assumption
 // digit is not local maximum 
            for (j=0; j<nn; j++) {
-               if (HitMap(cath)->TestHit(x[j], y[j])==kEmpty) continue;
-               digt=(AliMUONDigit*) HitMap(cath)->GetHit(x[j], y[j]);
-               isec=Segmentation(cath)->Sector(x[j], y[j]);
-               Float_t a1 = Segmentation(cath)->Dpx(isec)*Segmentation(cath)->Dpy(isec);
+               if (fHitMap[cath]->TestHit(x[j], y[j])==kEmpty) continue;
+               digt=(AliMUONDigit*) fHitMap[cath]->GetHit(x[j], y[j]);
+               isec=fInput->Segmentation(cath)->Sector(x[j], y[j]);
+               Float_t a1 = fInput->Segmentation(cath)->Dpx(isec)*fInput->Segmentation(cath)->Dpy(isec);
                if (digt->fSignal/a1 > fQ[i][cath]/a0) {
                    isLocal[i][cath]=kFALSE;
                    break;
@@ -938,9 +895,9 @@ void AliMUONClusterFinderVS::FindLocalMaxima(AliMUONRawCluster* c)
        cath1=1;
        
        for (i=0; i<fMul[cath]; i++) {
-           isec=Segmentation(cath)->Sector(fIx[i][cath],fIy[i][cath]);
-           dpy=Segmentation(cath)->Dpy(isec);
-           dpx=Segmentation(cath)->Dpx(isec);
+           isec=fInput->Segmentation(cath)->Sector(fIx[i][cath],fIy[i][cath]);
+           dpy=fInput->Segmentation(cath)->Dpy(isec);
+           dpx=fInput->Segmentation(cath)->Dpx(isec);
            if (isLocal[i][cath]) continue;
 // Pad position should be consistent with position of local maxima on the opposite cathode
            if ((TMath::Abs(fX[i][cath]-fX[fIndLocal[0][cath1]][cath1]) > dpx/2.) && 
@@ -950,18 +907,18 @@ void AliMUONClusterFinderVS::FindLocalMaxima(AliMUONRawCluster* c)
 // get neighbours for that digit and assume that it is local maximum       
            isLocal[i][cath]=kTRUE;
 // compare signal to that on the two neighbours on the left and on the right
-           Segmentation(cath)->GetPadIxy(fX[i][cath],fY[i][cath]+dpy,ix,iy);
+           fInput->Segmentation(cath)->GetPadI(fX[i][cath],fY[i][cath]+dpy,0,ix,iy);
 // iNN counts the number of neighbours with signal, it should be 1 or 2
            Int_t iNN=0;
-           if (HitMap(cath)->TestHit(ix, iy)!=kEmpty) {
+           if (fHitMap[cath]->TestHit(ix, iy)!=kEmpty) {
                iNN++;
-               digt=(AliMUONDigit*) HitMap(cath)->GetHit(ix,iy);
+               digt=(AliMUONDigit*) fHitMap[cath]->GetHit(ix,iy);
                if (digt->fSignal > fQ[i][cath]) isLocal[i][cath]=kFALSE;
            }
-           Segmentation(cath)->GetPadIxy(fX[i][cath],fY[i][cath]-dpy,ix,iy);
-           if (HitMap(cath)->TestHit(ix, iy)!=kEmpty) {
+           fInput->Segmentation(cath)->GetPadI(fX[i][cath],fY[i][cath]-dpy,0,ix,iy);
+           if (fHitMap[cath]->TestHit(ix, iy)!=kEmpty) {
                iNN++;
-               digt=(AliMUONDigit*) HitMap(cath)->GetHit(ix,iy);
+               digt=(AliMUONDigit*) fHitMap[cath]->GetHit(ix,iy);
                if (digt->fSignal > fQ[i][cath]) isLocal[i][cath]=kFALSE;
            }
            if (isLocal[i][cath] && iNN>0) {
@@ -991,9 +948,9 @@ void AliMUONClusterFinderVS::FindLocalMaxima(AliMUONRawCluster* c)
 //
 //  Loop over cluster digits
        for (i=0; i<fMul[cath]; i++) {
-           isec=Segmentation(cath)->Sector(fIx[i][cath],fIy[i][cath]);
-           dpx=Segmentation(cath)->Dpx(isec);
-           dpy=Segmentation(cath)->Dpy(isec);
+           isec=fInput->Segmentation(cath)->Sector(fIx[i][cath],fIy[i][cath]);
+           dpx=fInput->Segmentation(cath)->Dpx(isec);
+           dpy=fInput->Segmentation(cath)->Dpy(isec);
            if (isLocal[i][cath]) continue;
 // Pad position should be consistent with position of local maxima on the opposite cathode
            if ((TMath::Abs(fY[i][cath]-fY[fIndLocal[0][cath1]][cath1]) > dpy/2.) && 
@@ -1003,18 +960,18 @@ void AliMUONClusterFinderVS::FindLocalMaxima(AliMUONRawCluster* c)
 // get neighbours for that digit and assume that it is local maximum       
            isLocal[i][cath]=kTRUE;
 // compare signal to that on the two neighbours on the left and on the right
-           Segmentation(cath)->GetPadIxy(fX[i][cath]+dpx,fY[i][cath],ix,iy);
+           fInput->Segmentation(cath)->GetPadI(fX[i][cath]+dpx,fY[i][cath],0,ix,iy);
 // iNN counts the number of neighbours with signal, it should be 1 or 2
            Int_t iNN=0;
-           if (HitMap(cath)->TestHit(ix, iy)!=kEmpty) {
+           if (fHitMap[cath]->TestHit(ix, iy)!=kEmpty) {
                iNN++;
-               digt=(AliMUONDigit*) HitMap(cath)->GetHit(ix,iy);
+               digt=(AliMUONDigit*) fHitMap[cath]->GetHit(ix,iy);
                if (digt->fSignal > fQ[i][cath]) isLocal[i][cath]=kFALSE;
            }
-           Segmentation(cath)->GetPadIxy(fX[i][cath]-dpx,fY[i][cath],ix,iy);
-           if (HitMap(cath)->TestHit(ix, iy)!=kEmpty) {
+           fInput->Segmentation(cath)->GetPadI(fX[i][cath]-dpx,fY[i][cath],0,ix,iy);
+           if (fHitMap[cath]->TestHit(ix, iy)!=kEmpty) {
                iNN++;
-               digt=(AliMUONDigit*) HitMap(cath)->GetHit(ix,iy);
+               digt=(AliMUONDigit*) fHitMap[cath]->GetHit(ix,iy);
                if (digt->fSignal > fQ[i][cath]) isLocal[i][cath]=kFALSE;
            }
            if (isLocal[i][cath] && iNN>0) {
@@ -1043,7 +1000,7 @@ void  AliMUONClusterFinderVS::FillCluster(AliMUONRawCluster* c, Int_t flag, Int_
 //  Completes cluster information starting from list of digits
 //
     AliMUONDigit* dig;
-    Float_t x, y;
+    Float_t x, y, z;
     Int_t  ix, iy;
     
     if (cath==1) {
@@ -1062,7 +1019,7 @@ void  AliMUONClusterFinderVS::FillCluster(AliMUONRawCluster* c, Int_t flag, Int_
 //    fprintf(stderr,"\n fPeakSignal %d\n",c->fPeakSignal[cath]);
     for (Int_t i=0; i<c->fMultiplicity[cath]; i++)
     {
-       dig= (AliMUONDigit*)Digits(cath)->UncheckedAt(c->fIndexMap[i][cath]);
+       dig= fInput->Digit(cath,c->fIndexMap[i][cath]);
        ix=dig->fPadX+c->fOffsetMap[i][cath];
        iy=dig->fPadY;
        Int_t q=dig->fSignal;
@@ -1086,7 +1043,7 @@ void  AliMUONClusterFinderVS::FillCluster(AliMUONRawCluster* c, Int_t flag, Int_
        }
 //
        if (flag) {
-           Segmentation(cath)->GetPadCxy(ix, iy, x, y);
+           fInput->Segmentation(cath)->GetPadC(ix, iy, x, y, z);
            c->fX[cath] += q*x;
            c->fY[cath] += q*y;
            c->fQ[cath] += q;
@@ -1097,20 +1054,20 @@ void  AliMUONClusterFinderVS::FillCluster(AliMUONRawCluster* c, Int_t flag, Int_
 
     if (flag) {
        c->fX[cath]/=c->fQ[cath];
-       c->fX[cath]=Segmentation(cath)->GetAnod(c->fX[cath]);
+       c->fX[cath]=fInput->Segmentation(cath)->GetAnod(c->fX[cath]);
        c->fY[cath]/=c->fQ[cath]; 
 //
 //  apply correction to the coordinate along the anode wire
 //
        x=c->fX[cath];   
        y=c->fY[cath];
-       Segmentation(cath)->GetPadIxy(x, y, ix, iy);
-       Segmentation(cath)->GetPadCxy(ix, iy, x, y);
-       Int_t isec=Segmentation(cath)->Sector(ix,iy);
-       TF1* cogCorr = Segmentation(cath)->CorrFunc(isec-1);
+       fInput->Segmentation(cath)->GetPadI(x, y, 0, ix, iy);
+       fInput->Segmentation(cath)->GetPadC(ix, iy, x, y, z);
+       Int_t isec=fInput->Segmentation(cath)->Sector(ix,iy);
+       TF1* cogCorr = fInput->Segmentation(cath)->CorrFunc(isec-1);
        
        if (cogCorr) {
-           Float_t yOnPad=(c->fY[cath]-y)/Segmentation(cath)->Dpy(isec);
+           Float_t yOnPad=(c->fY[cath]-y)/fInput->Segmentation(cath)->Dpy(isec);
            c->fY[cath]=c->fY[cath]-cogCorr->Eval(yOnPad, 0, 0);
        }
     }
@@ -1129,14 +1086,14 @@ void  AliMUONClusterFinderVS::FillCluster(AliMUONRawCluster* c, Int_t cath)
        dr0 = 10000;
     }
     
-    Float_t xpad, ypad;
+    Float_t xpad, ypad, zpad;
     Float_t dx, dy, dr;
 
     for (Int_t i=0; i<c->fMultiplicity[cath]; i++)
     {
-      dig= (AliMUONDigit*)Digits(cath)->UncheckedAt(c->fIndexMap[i][cath]);
-      Segmentation(cath)->
-       GetPadCxy(dig->fPadX,dig->fPadY,xpad,ypad);
+       dig = fInput->Digit(cath,c->fIndexMap[i][cath]);
+       fInput->Segmentation(cath)->
+       GetPadC(dig->fPadX,dig->fPadY,xpad,ypad, zpad);
        fprintf(stderr,"x %f y %f cx %f cy %f\n",xpad,ypad,c->fX[0],c->fY[0]);
        dx = xpad - c->fX[0];
        dy = ypad - c->fY[0];
@@ -1161,7 +1118,7 @@ void  AliMUONClusterFinderVS::FillCluster(AliMUONRawCluster* c, Int_t cath)
     } // loop over digits
 
 //  apply correction to the coordinate along the anode wire
-    c->fX[cath]=Segmentation(cath)->GetAnod(c->fX[cath]);
+    c->fX[cath]=fInput->Segmentation(cath)->GetAnod(c->fX[cath]);
 }
 
 void  AliMUONClusterFinderVS::FindCluster(Int_t i, Int_t j, Int_t cath, AliMUONRawCluster &c){
@@ -1172,8 +1129,8 @@ void  AliMUONClusterFinderVS::FindCluster(Int_t i, Int_t j, Int_t cath, AliMUONR
 //  Add i,j as element of the cluster
 //
 
-    Int_t idx = HitMap(cath)->GetHitIndex(i,j);
-    AliMUONDigit* dig = (AliMUONDigit*) HitMap(cath)->GetHit(i,j);
+    Int_t idx = fHitMap[cath]->GetHitIndex(i,j);
+    AliMUONDigit* dig = (AliMUONDigit*) fHitMap[cath]->GetHit(i,j);
     Int_t q=dig->fSignal;
     Int_t theX=dig->fPadX;
     Int_t theY=dig->fPadY;    
@@ -1198,12 +1155,9 @@ void  AliMUONClusterFinderVS::FindCluster(Int_t i, Int_t j, Int_t cath, AliMUONR
     if (mu > 0) {
        for (Int_t ind=mu-1; ind>=0; ind--) {
            Int_t ist=(c.fIndexMap)[ind][cath];
-           Int_t ql=((AliMUONDigit*)Digits(cath)
-                     ->UncheckedAt(ist))->fSignal;
-           Int_t ix=((AliMUONDigit*)Digits(cath)
-                     ->UncheckedAt(ist))->fPadX;
-           Int_t iy=((AliMUONDigit*)Digits(cath)
-                     ->UncheckedAt(ist))->fPadY;
+           Int_t ql=fInput->Digit(cath, ist)->fSignal;
+           Int_t ix=fInput->Digit(cath, ist)->fPadX;
+           Int_t iy=fInput->Digit(cath, ist)->fPadY;
 
            if (q>ql || (q==ql && theX > ix && theY < iy)) {
                c.fIndexMap[ind][cath]=idx;
@@ -1221,53 +1175,53 @@ void  AliMUONClusterFinderVS::FindCluster(Int_t i, Int_t j, Int_t cath, AliMUONR
     }
 
 // Prepare center of gravity calculation
-    Float_t x, y;
-    Segmentation(cath)->GetPadCxy(i, j, x, y);
+    Float_t x, y, z;
+    fInput->Segmentation(cath)->GetPadC(i, j, x, y, z);
            
     c.fX[cath] += q*x;
     c.fY[cath] += q*y;
     c.fQ[cath] += q;
 // Flag hit as taken  
-    HitMap(cath)->FlagHit(i,j);
+    fHitMap[cath]->FlagHit(i,j);
 //
 //  Now look recursively for all neighbours and pad hit on opposite cathode
 //
 //  Loop over neighbours
     Int_t ix,iy;
     Int_t nn;
-    Int_t xList[kMaxNeighbours], yList[kMaxNeighbours];
-    Segmentation(cath)->Neighbours(i,j,&nn,xList,yList);
+    Int_t xList[10], yList[10];
+    fInput->Segmentation(cath)->Neighbours(i,j,&nn,xList,yList);
     for (Int_t in=0; in<nn; in++) {
        ix=xList[in];
        iy=yList[in];
-       if (HitMap(cath)->TestHit(ix,iy)==kUnused) FindCluster(ix, iy, cath, c);
+       if (fHitMap[cath]->TestHit(ix,iy)==kUnused) FindCluster(ix, iy, cath, c);
     }
 //  Neighbours on opposite cathode 
 //  Take into account that several pads can overlap with the present pad
     Float_t xmin, xmax, ymin, ymax, xc, yc;
     Int_t iop;
-    Int_t isec=Segmentation(cath)->Sector(i,j);    
+    Int_t isec=fInput->Segmentation(cath)->Sector(i,j);    
     if (cath==0) {
        iop=1;
-       xmin=x-Segmentation(cath)->Dpx(isec);
-       xmax=x+Segmentation(cath)->Dpx(isec);
+       xmin=x-fInput->Segmentation(cath)->Dpx(isec);
+       xmax=x+fInput->Segmentation(cath)->Dpx(isec);
        xc=xmin+.001;
        while (xc < xmax) {
-           xc+=Segmentation(iop)->Dpx(isec);
-           Segmentation(iop)->GetPadIxy(xc,y,ix,iy);
-           if (ix>=(Segmentation(iop)->Npx()) || (iy>=Segmentation(iop)->Npy())) continue;
-           if (HitMap(iop)->TestHit(ix,iy)==kUnused) FindCluster(ix, iy, iop, c);
+           xc+=fInput->Segmentation(iop)->Dpx(isec);
+           fInput->Segmentation(iop)->GetPadI(xc,y,0,ix,iy);
+           if (ix>=(fInput->Segmentation(iop)->Npx()) || (iy>=fInput->Segmentation(iop)->Npy())) continue;
+           if (fHitMap[iop]->TestHit(ix,iy)==kUnused) FindCluster(ix, iy, iop, c);
        }
     } else {
        iop=0;
-       ymin=y-Segmentation(cath)->Dpy(isec);
-       ymax=y+Segmentation(cath)->Dpy(isec);
+       ymin=y-fInput->Segmentation(cath)->Dpy(isec);
+       ymax=y+fInput->Segmentation(cath)->Dpy(isec);
        yc=ymin+.001;
        while (yc < ymax) {
-           yc+=Segmentation(iop)->Dpy(isec);
-           Segmentation(iop)->GetPadIxy(x,yc,ix,iy);
-           if (ix>=(Segmentation(iop)->Npx()) || (iy>=Segmentation(iop)->Npy())) continue;
-           if (HitMap(iop)->TestHit(ix,iy)==kUnused) FindCluster(ix, iy, iop, c);
+           yc+=fInput->Segmentation(iop)->Dpy(isec);
+           fInput->Segmentation(iop)->GetPadI(x,yc,0,ix,iy);
+           if (ix>=(fInput->Segmentation(iop)->Npx()) || (iy>=fInput->Segmentation(iop)->Npy())) continue;
+           if (fHitMap[iop]->TestHit(ix,iy)==kUnused) FindCluster(ix, iy, iop, c);
        }
     }
 }
@@ -1281,26 +1235,26 @@ void AliMUONClusterFinderVS::FindRawClusters()
   // fills the tree with raw clusters
   //
 
-    if (!NDigits(0) && !NDigits(1)) return;
+    if (!fInput->NDigits(0) && !fInput->NDigits(1)) return;
 
-    fHitMap  = new AliMUONHitMapA1(fSegmentation , fDigits);
-    fHitMap2 = new AliMUONHitMapA1(fSegmentation2, fDigits2);
+    fHitMap[0]  = new AliMUONHitMapA1(fInput->Segmentation(0), fInput->Digits(0));
+    fHitMap[1]  = new AliMUONHitMapA1(fInput->Segmentation(1), fInput->Digits(1));
 
     AliMUONDigit *dig;
 
     Int_t ndig, cath;
     Int_t nskip=0;
     Int_t ncls=0;
-    HitMap(0)->FillHits();
-    HitMap(1)->FillHits();
+    fHitMap[0]->FillHits();
+    fHitMap[1]->FillHits();
 //
 //  Outer Loop over Cathodes
     for (cath=0; cath<2; cath++) {
-       for (ndig=0; ndig<NDigits(cath); ndig++) {
-           dig = (AliMUONDigit*)Digits(cath)->UncheckedAt(ndig);
+       for (ndig=0; ndig<fInput->NDigits(cath); ndig++) {
+           dig = fInput->Digit(cath, ndig);
            Int_t i=dig->fPadX;
            Int_t j=dig->fPadY;
-           if (HitMap(cath)->TestHit(i,j)==kUsed ||fHitMap->TestHit(i,j)==kEmpty) {
+           if (fHitMap[cath]->TestHit(i,j)==kUsed ||fHitMap[0]->TestHit(i,j)==kEmpty) {
                nskip++;
                continue;
            }
@@ -1319,10 +1273,10 @@ void AliMUONClusterFinderVS::FindRawClusters()
 
            // center of gravity
            c.fX[0] /= c.fQ[0];
-           c.fX[0]=Segmentation(0)->GetAnod(c.fX[0]);
+           c.fX[0]=fInput->Segmentation(0)->GetAnod(c.fX[0]);
            c.fY[0] /= c.fQ[0];
            c.fX[1] /= c.fQ[1];
-           c.fX[1]=Segmentation(0)->GetAnod(c.fX[1]);
+           c.fX[1]=fInput->Segmentation(0)->GetAnod(c.fX[1]);
            c.fY[1] /= c.fQ[1];
            fprintf(stderr,"\n Cathode 1 multiplicite %d X(CG) %f Y(CG) %f\n",c.fMultiplicity[0],c.fX[0],c.fY[0]);
            fprintf(stderr," Cathode 2 multiplicite %d X(CG) %f Y(CG) %f\n",c.fMultiplicity[1],c.fX[1],c.fY[1]);
@@ -1332,9 +1286,9 @@ void AliMUONClusterFinderVS::FindRawClusters()
            Bool_t fitted;
            
            fitted=SingleMathiesonFit(&c, 0);
-           c.fX[0]=Segmentation(0)->GetAnod(c.fX[0]);
+           c.fX[0]=fInput->Segmentation(0)->GetAnod(c.fX[0]);
            fitted=SingleMathiesonFit(&c, 1);       
-           c.fX[1]=Segmentation(1)->GetAnod(c.fX[1]);
+           c.fX[1]=fInput->Segmentation(1)->GetAnod(c.fX[1]);
 */ 
 //
 //      Analyse cluster and decluster if necessary
@@ -1364,51 +1318,21 @@ void AliMUONClusterFinderVS::FindRawClusters()
        
        } // end loop ndig
     } // end loop cathodes
-    delete fHitMap;
-    delete fHitMap2;
+    delete fHitMap[0];
+    delete fHitMap[1];
 }
 
 Float_t AliMUONClusterFinderVS::SingleMathiesonFit(AliMUONRawCluster *c, Int_t cath)
 {
 //
-//  Initialise global variables for fit
-    Int_t i;
-    fMul[cath]=c->fMultiplicity[cath];
-    fgSegmentation[0]=Segmentation(cath);
-    fgResponse    =fResponse;
-    fgNbins[0]=fMul[cath];
-    Float_t qtot=0;
-//
-//  dump digit information into arrays
-//
-    for (i=0; i<fMul[cath]; i++)
-    {
-       fDig[i][cath]= (AliMUONDigit*)Digits(cath)->UncheckedAt(c->fIndexMap[i][cath]);
-       fIx[i][cath]= fDig[i][cath]->fPadX;
-       fIy[i][cath]= fDig[i][cath]->fPadY;
-       fQ[i][cath] = fDig[i][cath]->fSignal;
-       Segmentation(cath)->GetPadCxy(fIx[i][cath], fIy[i][cath], fX[i][cath], fY[i][cath]);
-       fgix[i][0]=fIx[i][cath];
-       fgiy[i][0]=fIy[i][cath];
-       fgCharge[i][0]=Float_t(fQ[i][cath]);
-       qtot+=fgCharge[i][0];
-    }
-
-    fgQtot[0]=qtot;
-    fgChargeTot[0]=Int_t(qtot);
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
     
-//
-    if (fgFirst) {
-       fgFirst=kFALSE;
-       fgMyMinuit = new TMinuit(5);
-    }
-
-    fgMyMinuit->SetFCN(fcnS1);
-    fgMyMinuit->mninit(2,10,7);
+    clusterInput.Fitter()->SetFCN(fcnS1);
+    clusterInput.Fitter()->mninit(2,10,7);
     Double_t arglist[20];
     Int_t ierflag=0;
     arglist[0]=1;
-//     fgMyMinuit->mnexcm("SET ERR",arglist,1,ierflag);
+//     clusterInput.Fitter()->mnexcm("SET ERR",arglist,1,ierflag);
 // Set starting values 
     static Double_t vstart[2];
     vstart[0]=c->fX[1];
@@ -1418,32 +1342,32 @@ Float_t AliMUONClusterFinderVS::SingleMathiesonFit(AliMUONRawCluster *c, Int_t c
 // lower and upper limits
     static Double_t lower[2], upper[2];
     Int_t ix,iy;
-    Segmentation(cath)->GetPadIxy(c->fX[cath], c->fY[cath], ix, iy);
-    Int_t isec=Segmentation(cath)->Sector(ix, iy);
-    lower[0]=vstart[0]-Segmentation(cath)->Dpx(isec)/2;
-    lower[1]=vstart[1]-Segmentation(cath)->Dpy(isec)/2;
+    fInput->Segmentation(cath)->GetPadI(c->fX[cath], c->fY[cath], 0, ix, iy);
+    Int_t isec=fInput->Segmentation(cath)->Sector(ix, iy);
+    lower[0]=vstart[0]-fInput->Segmentation(cath)->Dpx(isec)/2;
+    lower[1]=vstart[1]-fInput->Segmentation(cath)->Dpy(isec)/2;
     
-    upper[0]=lower[0]+Segmentation(cath)->Dpx(isec);
-    upper[1]=lower[1]+Segmentation(cath)->Dpy(isec);
+    upper[0]=lower[0]+fInput->Segmentation(cath)->Dpx(isec);
+    upper[1]=lower[1]+fInput->Segmentation(cath)->Dpy(isec);
     
 // step sizes
     static Double_t step[2]={0.0005, 0.0005};
     
-    fgMyMinuit->mnparm(0,"x1",vstart[0],step[0],lower[0],upper[0],ierflag);
-    fgMyMinuit->mnparm(1,"y1",vstart[1],step[1],lower[1],upper[1],ierflag);
+    clusterInput.Fitter()->mnparm(0,"x1",vstart[0],step[0],lower[0],upper[0],ierflag);
+    clusterInput.Fitter()->mnparm(1,"y1",vstart[1],step[1],lower[1],upper[1],ierflag);
 // ready for minimisation      
-    fgMyMinuit->SetPrintLevel(1);
-    fgMyMinuit->mnexcm("SET OUT", arglist, 0, ierflag);
+    clusterInput.Fitter()->SetPrintLevel(1);
+    clusterInput.Fitter()->mnexcm("SET OUT", arglist, 0, ierflag);
     arglist[0]= -1;
     arglist[1]= 0;
     
-    fgMyMinuit->mnexcm("SET NOGR", arglist, 0, ierflag);
-    fgMyMinuit->mnexcm("MIGRAD", arglist, 0, ierflag);
-    fgMyMinuit->mnexcm("EXIT" , arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("SET NOGR", arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("MIGRAD", arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("EXIT" , arglist, 0, ierflag);
     Double_t fmin, fedm, errdef;
     Int_t   npari, nparx, istat;
       
-    fgMyMinuit->mnstat(fmin, fedm, errdef, npari, nparx, istat);  
+    clusterInput.Fitter()->mnstat(fmin, fedm, errdef, npari, nparx, istat);  
     fFitStat=istat;
     
 // Print results
@@ -1452,8 +1376,8 @@ Float_t AliMUONClusterFinderVS::SingleMathiesonFit(AliMUONRawCluster *c, Int_t c
     TString chname;
     Double_t epxz, b1, b2;
     Int_t ierflg;
-    fgMyMinuit->mnpout(0, chname, xrec, epxz, b1, b2, ierflg); 
-    fgMyMinuit->mnpout(1, chname, yrec, epxz, b1, b2, ierflg); 
+    clusterInput.Fitter()->mnpout(0, chname, xrec, epxz, b1, b2, ierflg);      
+    clusterInput.Fitter()->mnpout(1, chname, yrec, epxz, b1, b2, ierflg);      
     fXFit[cath]=xrec;
     fYFit[cath]=yrec;
     return fmin;
@@ -1463,13 +1387,9 @@ Float_t AliMUONClusterFinderVS::CombiSingleMathiesonFit(AliMUONRawCluster *c)
 {
 // Perform combined Mathieson fit on both cathode planes
 //
-    if (fgFirst) {
-       fgFirst=kFALSE;
-       fgMyMinuit = new TMinuit(5);
-    }
-
-    fgMyMinuit->SetFCN(fcnCombiS1);
-    fgMyMinuit->mninit(2,10,7);
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
+    clusterInput.Fitter()->SetFCN(fcnCombiS1);
+    clusterInput.Fitter()->mninit(2,10,7);
     Double_t arglist[20];
     Int_t ierflag=0;
     arglist[0]=1;
@@ -1481,12 +1401,12 @@ Float_t AliMUONClusterFinderVS::CombiSingleMathiesonFit(AliMUONRawCluster *c)
 // lower and upper limits
     static Double_t lower[2], upper[2];
     Int_t ix,iy,isec;
-    Segmentation(0)->GetPadIxy(fXInit[0], fYInit[0], ix, iy);
-    isec=Segmentation(0)->Sector(ix, iy);
-    Float_t dpy=Segmentation(0)->Dpy(isec)/2;
-    Segmentation(1)->GetPadIxy(fXInit[0], fYInit[0], ix, iy);
-    isec=Segmentation(1)->Sector(ix, iy);
-    Float_t dpx=Segmentation(1)->Dpx(isec)/2;
+    fInput->Segmentation(0)->GetPadI(fXInit[0], fYInit[0], 0, ix, iy);
+    isec=fInput->Segmentation(0)->Sector(ix, iy);
+    Float_t dpy=fInput->Segmentation(0)->Dpy(isec)/2;
+    fInput->Segmentation(1)->GetPadI(fXInit[0], fYInit[0], 0, ix, iy);
+    isec=fInput->Segmentation(1)->Sector(ix, iy);
+    Float_t dpx=fInput->Segmentation(1)->Dpx(isec)/2;
 
 
     lower[0]=vstart[0]-dpx;
@@ -1498,21 +1418,21 @@ Float_t AliMUONClusterFinderVS::CombiSingleMathiesonFit(AliMUONRawCluster *c)
 // step sizes
     static Double_t step[2]={0.00001, 0.0001};
     
-    fgMyMinuit->mnparm(0,"x1",vstart[0],step[0],lower[0],upper[0],ierflag);
-    fgMyMinuit->mnparm(1,"y1",vstart[1],step[1],lower[1],upper[1],ierflag);
+    clusterInput.Fitter()->mnparm(0,"x1",vstart[0],step[0],lower[0],upper[0],ierflag);
+    clusterInput.Fitter()->mnparm(1,"y1",vstart[1],step[1],lower[1],upper[1],ierflag);
 // ready for minimisation      
-    fgMyMinuit->SetPrintLevel(1);
-    fgMyMinuit->mnexcm("SET OUT", arglist, 0, ierflag);
+    clusterInput.Fitter()->SetPrintLevel(1);
+    clusterInput.Fitter()->mnexcm("SET OUT", arglist, 0, ierflag);
     arglist[0]= -1;
     arglist[1]= 0;
     
-    fgMyMinuit->mnexcm("SET NOGR", arglist, 0, ierflag);
-    fgMyMinuit->mnexcm("MIGRAD", arglist, 0, ierflag);
-    fgMyMinuit->mnexcm("EXIT" , arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("SET NOGR", arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("MIGRAD", arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("EXIT" , arglist, 0, ierflag);
     Double_t fmin, fedm, errdef;
     Int_t   npari, nparx, istat;
       
-    fgMyMinuit->mnstat(fmin, fedm, errdef, npari, nparx, istat);  
+    clusterInput.Fitter()->mnstat(fmin, fedm, errdef, npari, nparx, istat);  
     fFitStat=istat;
     
 // Print results
@@ -1521,8 +1441,8 @@ Float_t AliMUONClusterFinderVS::CombiSingleMathiesonFit(AliMUONRawCluster *c)
     TString chname;
     Double_t epxz, b1, b2;
     Int_t ierflg;
-    fgMyMinuit->mnpout(0, chname, xrec, epxz, b1, b2, ierflg); 
-    fgMyMinuit->mnpout(1, chname, yrec, epxz, b1, b2, ierflg); 
+    clusterInput.Fitter()->mnpout(0, chname, xrec, epxz, b1, b2, ierflg);      
+    clusterInput.Fitter()->mnpout(1, chname, yrec, epxz, b1, b2, ierflg);      
     fXFit[0]=xrec;
     fYFit[0]=yrec;
     return fmin;
@@ -1532,29 +1452,9 @@ Bool_t AliMUONClusterFinderVS::DoubleMathiesonFit(AliMUONRawCluster *c, Int_t ca
 {
 //
 //  Initialise global variables for fit
-    Int_t i,j;
-    
-    fgSegmentation[0]=Segmentation(cath);
-    fgResponse    =fResponse;
-    fgNbins[0]=fMul[cath];
-    Float_t qtot=0;
-       
-    for (i=0; i<fMul[cath]; i++) {
-       fgix[i][0]=fIx[i][cath];
-       fgiy[i][0]=fIy[i][cath];
-       fgCharge[i][0]=Float_t(fQ[i][cath]);
-       qtot+=fgCharge[i][0];
-    }
-    fgQtot[0]=qtot;
-    fgChargeTot[0]=Int_t(qtot);
-    
-//
-    if (fgFirst) {
-       fgFirst=kFALSE;
-       fgMyMinuit = new TMinuit(5);
-    }
-    fgMyMinuit->SetFCN(fcnS2);
-    fgMyMinuit->mninit(5,10,7);
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
+    clusterInput.Fitter()->SetFCN(fcnS2);
+    clusterInput.Fitter()->mninit(5,10,7);
     Double_t arglist[20];
     Int_t ierflag=0;
     arglist[0]=1;
@@ -1568,92 +1468,55 @@ Bool_t AliMUONClusterFinderVS::DoubleMathiesonFit(AliMUONRawCluster *c, Int_t ca
        Float_t(fQ[fIndLocal[0][cath]][cath]+fQ[fIndLocal[1][cath]][cath]);
 // lower and upper limits
     static Double_t lower[5], upper[5];
-    Int_t isec=Segmentation(cath)->Sector(fIx[fIndLocal[0][cath]][cath], fIy[fIndLocal[0][cath]][cath]);
-    lower[0]=vstart[0]-Segmentation(cath)->Dpx(isec);
-    lower[1]=vstart[1]-Segmentation(cath)->Dpy(isec);
+    Int_t isec=fInput->Segmentation(cath)->Sector(fIx[fIndLocal[0][cath]][cath], fIy[fIndLocal[0][cath]][cath]);
+    lower[0]=vstart[0]-fInput->Segmentation(cath)->Dpx(isec);
+    lower[1]=vstart[1]-fInput->Segmentation(cath)->Dpy(isec);
     
-    upper[0]=lower[0]+2.*Segmentation(cath)->Dpx(isec);
-    upper[1]=lower[1]+2.*Segmentation(cath)->Dpy(isec);
+    upper[0]=lower[0]+2.*fInput->Segmentation(cath)->Dpx(isec);
+    upper[1]=lower[1]+2.*fInput->Segmentation(cath)->Dpy(isec);
     
-    isec=Segmentation(cath)->Sector(fIx[fIndLocal[1][cath]][cath], fIy[fIndLocal[1][cath]][cath]);
-    lower[2]=vstart[2]-Segmentation(cath)->Dpx(isec)/2;
-    lower[3]=vstart[3]-Segmentation(cath)->Dpy(isec)/2;
+    isec=fInput->Segmentation(cath)->Sector(fIx[fIndLocal[1][cath]][cath], fIy[fIndLocal[1][cath]][cath]);
+    lower[2]=vstart[2]-fInput->Segmentation(cath)->Dpx(isec)/2;
+    lower[3]=vstart[3]-fInput->Segmentation(cath)->Dpy(isec)/2;
     
-    upper[2]=lower[2]+Segmentation(cath)->Dpx(isec);
-    upper[3]=lower[3]+Segmentation(cath)->Dpy(isec);
+    upper[2]=lower[2]+fInput->Segmentation(cath)->Dpx(isec);
+    upper[3]=lower[3]+fInput->Segmentation(cath)->Dpy(isec);
     
     lower[4]=0.;
     upper[4]=1.;
 // step sizes
     static Double_t step[5]={0.0005, 0.0005, 0.0005, 0.0005, 0.0001};
     
-    fgMyMinuit->mnparm(0,"x1",vstart[0],step[0],lower[0],upper[0],ierflag);
-    fgMyMinuit->mnparm(1,"y1",vstart[1],step[1],lower[1],upper[1],ierflag);
-    fgMyMinuit->mnparm(2,"x2",vstart[2],step[2],lower[2],upper[2],ierflag);
-    fgMyMinuit->mnparm(3,"y2",vstart[3],step[3],lower[3],upper[3],ierflag);
-    fgMyMinuit->mnparm(4,"a0",vstart[4],step[4],lower[4],upper[4],ierflag);
+    clusterInput.Fitter()->mnparm(0,"x1",vstart[0],step[0],lower[0],upper[0],ierflag);
+    clusterInput.Fitter()->mnparm(1,"y1",vstart[1],step[1],lower[1],upper[1],ierflag);
+    clusterInput.Fitter()->mnparm(2,"x2",vstart[2],step[2],lower[2],upper[2],ierflag);
+    clusterInput.Fitter()->mnparm(3,"y2",vstart[3],step[3],lower[3],upper[3],ierflag);
+    clusterInput.Fitter()->mnparm(4,"a0",vstart[4],step[4],lower[4],upper[4],ierflag);
 // ready for minimisation      
-    fgMyMinuit->SetPrintLevel(-1);
-    fgMyMinuit->mnexcm("SET OUT", arglist, 0, ierflag);
+    clusterInput.Fitter()->SetPrintLevel(-1);
+    clusterInput.Fitter()->mnexcm("SET OUT", arglist, 0, ierflag);
     arglist[0]= -1;
     arglist[1]= 0;
     
-    fgMyMinuit->mnexcm("SET NOGR", arglist, 0, ierflag);
-    fgMyMinuit->mnexcm("MIGRAD", arglist, 0, ierflag);
-    fgMyMinuit->mnexcm("EXIT" , arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("SET NOGR", arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("MIGRAD", arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("EXIT" , arglist, 0, ierflag);
 // Get fitted parameters
     Double_t xrec[2], yrec[2], qfrac;
     TString chname;
     Double_t epxz, b1, b2;
     Int_t ierflg;
-    fgMyMinuit->mnpout(0, chname, xrec[0], epxz, b1, b2, ierflg);      
-    fgMyMinuit->mnpout(1, chname, yrec[0], epxz, b1, b2, ierflg);      
-    fgMyMinuit->mnpout(2, chname, xrec[1], epxz, b1, b2, ierflg);      
-    fgMyMinuit->mnpout(3, chname, yrec[1], epxz, b1, b2, ierflg);      
-    fgMyMinuit->mnpout(4, chname, qfrac,   epxz, b1, b2, ierflg);      
+    clusterInput.Fitter()->mnpout(0, chname, xrec[0], epxz, b1, b2, ierflg);   
+    clusterInput.Fitter()->mnpout(1, chname, yrec[0], epxz, b1, b2, ierflg);   
+    clusterInput.Fitter()->mnpout(2, chname, xrec[1], epxz, b1, b2, ierflg);   
+    clusterInput.Fitter()->mnpout(3, chname, yrec[1], epxz, b1, b2, ierflg);   
+    clusterInput.Fitter()->mnpout(4, chname, qfrac,   epxz, b1, b2, ierflg);   
 
     Double_t fmin, fedm, errdef;
     Int_t   npari, nparx, istat;
       
-    fgMyMinuit->mnstat(fmin, fedm, errdef, npari, nparx, istat);  
+    clusterInput.Fitter()->mnstat(fmin, fedm, errdef, npari, nparx, istat);  
     fFitStat=istat;
-    
-    
-//
-// One cluster for each maximum
-//
-    for (j=0; j<2; j++) {
-       AliMUONRawCluster cnew;
-       cnew.fChi2[0]=Float_t(fmin);
-       
-       if (fNPeaks == 0) {
-           cnew.fNcluster[0]=-1;
-           cnew.fNcluster[1]=fNRawClusters;
-       } else {
-           cnew.fNcluster[0]=fNPeaks;
-           cnew.fNcluster[1]=0;
-       }
-       cnew.fMultiplicity[0]=0;
-       cnew.fX[0]=Float_t(xrec[j]);
-       cnew.fY[0]=Float_t(yrec[j]);
-       if (j==0) {
-           cnew.fQ[0]=Int_t(fgChargeTot[0]*qfrac);
-       } else {
-           cnew.fQ[0]=Int_t(fgChargeTot[0]*(1-qfrac));
-       }
-       fgSegmentation[0]->SetHit(xrec[j],yrec[j]);
-       for (i=0; i<fMul[cath]; i++) {
-           cnew.fIndexMap[cnew.fMultiplicity[0]][cath]=c->fIndexMap[i][cath];
-           fgSegmentation[0]->SetPad(fgix[i][0], fgiy[i][0]);
-           Float_t q1=fgResponse->IntXY(fgSegmentation[0]);
-           cnew.fContMap[cnew.fMultiplicity[0]][0]=(q1*cnew.fQ[0])/Float_t(fQ[i][cath]);
-           cnew.fMultiplicity[0]++;
-       }
-       FillCluster(&cnew,0,0);
-       cnew.fClusterType=cnew.PhysicsContribution();
-       AddRawCluster(cnew);
-       fNPeaks++;
-    }
     return kTRUE;
 }
 
@@ -1662,12 +1525,9 @@ Float_t AliMUONClusterFinderVS::CombiDoubleMathiesonFit(AliMUONRawCluster *c)
 //
 // Perform combined double Mathieson fit on both cathode planes
 //
-    if (fgFirst) {
-       fgFirst=kFALSE;
-       fgMyMinuit = new TMinuit(5);
-    }
-    fgMyMinuit->SetFCN(fcnCombiS2);
-    fgMyMinuit->mninit(6,10,7);
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
+    clusterInput.Fitter()->SetFCN(fcnCombiS2);
+    clusterInput.Fitter()->mninit(6,10,7);
     Double_t arglist[20];
     Int_t ierflag=0;
     arglist[0]=1;
@@ -1684,13 +1544,13 @@ Float_t AliMUONClusterFinderVS::CombiDoubleMathiesonFit(AliMUONRawCluster *c)
     Int_t ix,iy,isec;
     Float_t dpx, dpy;
     
-    Segmentation(1)->GetPadIxy(fXInit[0], fYInit[0], ix, iy);
-    isec=Segmentation(1)->Sector(ix, iy);
-    dpx=Segmentation(1)->Dpx(isec);
+    fInput->Segmentation(1)->GetPadI(fXInit[0], fYInit[0], 0, ix, iy);
+    isec=fInput->Segmentation(1)->Sector(ix, iy);
+    dpx=fInput->Segmentation(1)->Dpx(isec);
 
-    Segmentation(0)->GetPadIxy(fXInit[0], fYInit[0], ix, iy);
-    isec=Segmentation(0)->Sector(ix, iy);
-    dpy=Segmentation(0)->Dpy(isec);
+    fInput->Segmentation(0)->GetPadI(fXInit[0], fYInit[0], 0, ix, iy);
+    isec=fInput->Segmentation(0)->Sector(ix, iy);
+    dpy=fInput->Segmentation(0)->Dpy(isec);
 
     lower[0]=vstart[0]-dpx;
     lower[1]=vstart[1]-dpy;
@@ -1698,12 +1558,12 @@ Float_t AliMUONClusterFinderVS::CombiDoubleMathiesonFit(AliMUONRawCluster *c)
     upper[1]=vstart[1]+dpy;
 
 
-    Segmentation(1)->GetPadIxy(fXInit[1], fYInit[1], ix, iy);
-    isec=Segmentation(1)->Sector(ix, iy);
-    dpx=Segmentation(1)->Dpx(isec);
-    Segmentation(0)->GetPadIxy(fXInit[1], fYInit[1], ix, iy);
-    isec=Segmentation(0)->Sector(ix, iy);
-    dpy=Segmentation(0)->Dpy(isec);
+    fInput->Segmentation(1)->GetPadI(fXInit[1], fYInit[1], 0, ix, iy);
+    isec=fInput->Segmentation(1)->Sector(ix, iy);
+    dpx=fInput->Segmentation(1)->Dpx(isec);
+    fInput->Segmentation(0)->GetPadI(fXInit[1], fYInit[1], 0, ix, iy);
+    isec=fInput->Segmentation(0)->Sector(ix, iy);
+    dpy=fInput->Segmentation(0)->Dpy(isec);
 
     lower[2]=vstart[2]-dpx;
     lower[3]=vstart[3]-dpy;
@@ -1719,36 +1579,36 @@ Float_t AliMUONClusterFinderVS::CombiDoubleMathiesonFit(AliMUONRawCluster *c)
 // step sizes
     static Double_t step[6]={0.0005, 0.0005, 0.0005, 0.0005, 0.001, 0.001};
     
-    fgMyMinuit->mnparm(0,"x1",vstart[0],step[0],lower[0],upper[0],ierflag);
-    fgMyMinuit->mnparm(1,"y1",vstart[1],step[1],lower[1],upper[1],ierflag);
-    fgMyMinuit->mnparm(2,"x2",vstart[2],step[2],lower[2],upper[2],ierflag);
-    fgMyMinuit->mnparm(3,"y2",vstart[3],step[3],lower[3],upper[3],ierflag);
-    fgMyMinuit->mnparm(4,"a0",vstart[4],step[4],lower[4],upper[4],ierflag);
-    fgMyMinuit->mnparm(5,"a1",vstart[5],step[5],lower[5],upper[5],ierflag);
+    clusterInput.Fitter()->mnparm(0,"x1",vstart[0],step[0],lower[0],upper[0],ierflag);
+    clusterInput.Fitter()->mnparm(1,"y1",vstart[1],step[1],lower[1],upper[1],ierflag);
+    clusterInput.Fitter()->mnparm(2,"x2",vstart[2],step[2],lower[2],upper[2],ierflag);
+    clusterInput.Fitter()->mnparm(3,"y2",vstart[3],step[3],lower[3],upper[3],ierflag);
+    clusterInput.Fitter()->mnparm(4,"a0",vstart[4],step[4],lower[4],upper[4],ierflag);
+    clusterInput.Fitter()->mnparm(5,"a1",vstart[5],step[5],lower[5],upper[5],ierflag);
 // ready for minimisation      
-    fgMyMinuit->SetPrintLevel(-1);
-    fgMyMinuit->mnexcm("SET OUT", arglist, 0, ierflag);
+    clusterInput.Fitter()->SetPrintLevel(-1);
+    clusterInput.Fitter()->mnexcm("SET OUT", arglist, 0, ierflag);
     arglist[0]= -1;
     arglist[1]= 0;
     
-    fgMyMinuit->mnexcm("SET NOGR", arglist, 0, ierflag);
-    fgMyMinuit->mnexcm("MIGRAD", arglist, 0, ierflag);
-    fgMyMinuit->mnexcm("EXIT" , arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("SET NOGR", arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("MIGRAD", arglist, 0, ierflag);
+    clusterInput.Fitter()->mnexcm("EXIT" , arglist, 0, ierflag);
 // Get fitted parameters
     TString chname;
     Double_t epxz, b1, b2;
     Int_t ierflg;
-    fgMyMinuit->mnpout(0, chname, fXFit[0],  epxz, b1, b2, ierflg);    
-    fgMyMinuit->mnpout(1, chname, fYFit[0],  epxz, b1, b2, ierflg);    
-    fgMyMinuit->mnpout(2, chname, fXFit[1],  epxz, b1, b2, ierflg);    
-    fgMyMinuit->mnpout(3, chname, fYFit[1],  epxz, b1, b2, ierflg);    
-    fgMyMinuit->mnpout(4, chname, fQrFit[0], epxz, b1, b2, ierflg);    
-    fgMyMinuit->mnpout(5, chname, fQrFit[1], epxz, b1, b2, ierflg);    
+    clusterInput.Fitter()->mnpout(0, chname, fXFit[0],  epxz, b1, b2, ierflg); 
+    clusterInput.Fitter()->mnpout(1, chname, fYFit[0],  epxz, b1, b2, ierflg); 
+    clusterInput.Fitter()->mnpout(2, chname, fXFit[1],  epxz, b1, b2, ierflg); 
+    clusterInput.Fitter()->mnpout(3, chname, fYFit[1],  epxz, b1, b2, ierflg); 
+    clusterInput.Fitter()->mnpout(4, chname, fQrFit[0], epxz, b1, b2, ierflg); 
+    clusterInput.Fitter()->mnpout(5, chname, fQrFit[1], epxz, b1, b2, ierflg); 
 
     Double_t fmin, fedm, errdef;
     Int_t   npari, nparx, istat;
       
-    fgMyMinuit->mnstat(fmin, fedm, errdef, npari, nparx, istat);  
+    clusterInput.Fitter()->mnstat(fmin, fedm, errdef, npari, nparx, istat);  
     fFitStat=istat;
     
     fChi2[0]=fmin;
@@ -1762,7 +1622,7 @@ void AliMUONClusterFinderVS::Split(AliMUONRawCluster* c)
 // One cluster for each maximum
 //
     Int_t i, j, cath;
-    
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
     for (j=0; j<2; j++) {
        AliMUONRawCluster cnew;
        for (cath=0; cath<2; cath++) {
@@ -1779,16 +1639,16 @@ void AliMUONClusterFinderVS::Split(AliMUONRawCluster* c)
            cnew.fX[cath]=Float_t(fXFit[j]);
            cnew.fY[cath]=Float_t(fYFit[j]);
            if (j==0) {
-               cnew.fQ[cath]=Int_t(fgChargeTot[cath]*fQrFit[cath]);
+               cnew.fQ[cath]=Int_t(clusterInput.TotalCharge(cath)*fQrFit[cath]);
            } else {
-               cnew.fQ[cath]=Int_t(fgChargeTot[cath]*(1-fQrFit[cath]));
+               cnew.fQ[cath]=Int_t(clusterInput.TotalCharge(cath)*(1-fQrFit[cath]));
            }
-           fgSegmentation[cath]->SetHit(fXFit[j],fYFit[j]);
+           fInput->Segmentation(cath)->SetHit(fXFit[j],fYFit[j],0);
            for (i=0; i<fMul[cath]; i++) {
                cnew.fIndexMap[cnew.fMultiplicity[cath]][cath]=
                    c->fIndexMap[i][cath];
-               fgSegmentation[cath]->SetPad(fgix[i][cath], fgiy[i][cath]);
-               Float_t q1=fgResponse->IntXY(fgSegmentation[cath]);
+               fInput->Segmentation(cath)->SetPad(fIx[i][cath], fIy[i][cath]);
+               Float_t q1=fInput->Response()->IntXY(fInput->Segmentation(cath));
                cnew.fContMap[i][cath]
                    =(q1*Float_t(cnew.fQ[cath]))/Float_t(fQ[i][cath]);
                cnew.fMultiplicity[cath]++;
@@ -1804,97 +1664,21 @@ void AliMUONClusterFinderVS::Split(AliMUONRawCluster* c)
 }
 
 
-Float_t DiscrChargeS1(Int_t i,Double_t *par) 
-{
-// par[0]    x-position of cluster
-// par[1]    y-position of cluster
-
-   fgSegmentation[0]->SetPad(fgix[i][0], fgiy[i][0]);
-//  First Cluster
-   fgSegmentation[0]->SetHit(par[0],par[1]);
-   Float_t q1=fgResponse->IntXY(fgSegmentation[0]);
-    
-   Float_t value = fgQtot[0]*q1;
-   return value;
-}
-
-Float_t DiscrChargeCombiS1(Int_t i,Double_t *par, Int_t cath) 
-{
-// par[0]    x-position of cluster
-// par[1]    y-position of cluster
-
-   fgSegmentation[cath]->SetPad(fgix[i][cath], fgiy[i][cath]);
-//  First Cluster
-   fgSegmentation[cath]->SetHit(par[0],par[1]);
-   Float_t q1=fgResponse->IntXY(fgSegmentation[cath]);
-    
-   Float_t value = fgQtot[cath]*q1;
-   return value;
-}
-
-
-Float_t DiscrChargeS2(Int_t i,Double_t *par) 
-{
-// par[0]    x-position of first  cluster
-// par[1]    y-position of first  cluster
-// par[2]    x-position of second cluster
-// par[3]    y-position of second cluster
-// par[4]    charge fraction of first  cluster
-// 1-par[4]  charge fraction of second cluster
-
-   fgSegmentation[0]->SetPad(fgix[i][0], fgiy[i][0]);
-//  First Cluster
-   fgSegmentation[0]->SetHit(par[0],par[1]);
-   Float_t q1=fgResponse->IntXY(fgSegmentation[0]);
-    
-//  Second Cluster
-   fgSegmentation[0]->SetHit(par[2],par[3]);
-   Float_t q2=fgResponse->IntXY(fgSegmentation[0]);
-    
-   Float_t value = fgQtot[0]*(par[4]*q1+(1.-par[4])*q2);
-   return value;
-}
-
-Float_t DiscrChargeCombiS2(Int_t i,Double_t *par, Int_t cath) 
-{
-// par[0]    x-position of first  cluster
-// par[1]    y-position of first  cluster
-// par[2]    x-position of second cluster
-// par[3]    y-position of second cluster
-// par[4]    charge fraction of first  cluster
-// 1-par[4]  charge fraction of second cluster
-
-   fgSegmentation[cath]->SetPad(fgix[i][cath], fgiy[i][cath]);
-//  First Cluster
-   fgSegmentation[cath]->SetHit(par[0],par[1]);
-   Float_t q1=fgResponse->IntXY(fgSegmentation[cath]);
-    
-//  Second Cluster
-   fgSegmentation[cath]->SetHit(par[2],par[3]);
-   Float_t q2=fgResponse->IntXY(fgSegmentation[cath]);
-   Float_t value;
-   if (cath==0) {
-       value = fgQtot[0]*(par[4]*q1+(1.-par[4])*q2);
-   } else {
-       value = fgQtot[1]*(par[5]*q1+(1.-par[5])*q2);
-   }
-   return value;
-}
-
 //
 // Minimisation functions
 // Single Mathieson
 void fcnS1(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
 {
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());    
     Int_t i;
     Float_t delta;
     Float_t chisq=0;
     Float_t qcont=0;
     Float_t qtot=0;
-    
-    for (i=0; i<fgNbins[0]; i++) {
-       Float_t q0=fgCharge[i][0];
-       Float_t q1=DiscrChargeS1(i,par);
+
+    for (i=0; i<clusterInput.Nmul(0); i++) {
+       Float_t q0=clusterInput.Charge(i,0);
+       Float_t q1=clusterInput.DiscrChargeS1(i,par);
        delta=(q0-q1)/q0;
        chisq+=delta*delta;
        qcont+=q1;
@@ -1905,6 +1689,7 @@ void fcnS1(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
 
 void fcnCombiS1(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
 {
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());    
     Int_t i, cath;
     Float_t delta;
     Float_t chisq=0;
@@ -1914,35 +1699,35 @@ void fcnCombiS1(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t if
 
     for (cath=0; cath<2; cath++) {
 //     chisq=0;
-       for (i=0; i<fgNbins[cath]; i++) {
-           Float_t q0=fgCharge[i][cath];
-           Float_t q1=DiscrChargeCombiS1(i,par,cath);
+       for (i=0; i<clusterInput.Nmul(cath); i++) {
+           Float_t q0=clusterInput.Charge(i,cath);
+           Float_t q1=clusterInput.DiscrChargeCombiS1(i,par,cath);
            //      delta=(q0-q1);
            delta=(q0-q1)/q0;
            chisq+=delta*delta;
            qcont+=q1;
            qtot+=q0;
        }
-//     if (cath == 0) chi2temp=chisq/fgNbins[cath];
+//     if (cath == 0) chi2temp=chisq/clusterInput.Nbins[cath];
     }
-//    chisq = chisq/fgNbins[1]+chi2temp; 
-    
+//    chisq = chisq/clusterInput.Nbins[1]+chi2temp; 
     f=chisq;
 }
 
 // Double Mathieson
 void fcnS2(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
 {
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());    
     Int_t i;
     Float_t delta;
     Float_t chisq=0;
     Float_t qcont=0;
     Float_t qtot=0;
     
-    for (i=0; i<fgNbins[0]; i++) {
+    for (i=0; i<clusterInput.Nmul(0); i++) {
 
-       Float_t q0=fgCharge[i][0];
-       Float_t q1=DiscrChargeS2(i,par);
+       Float_t q0=clusterInput.Charge(i,0);
+       Float_t q1=clusterInput.DiscrChargeS2(i,par);
        delta=(q0-q1)/q0;
        chisq+=delta*delta;
        qcont+=q1;
@@ -1955,6 +1740,7 @@ void fcnS2(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
 // Double Mathieson
 void fcnCombiS2(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
 {
+    AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());    
     Int_t i, cath;
     Float_t delta;
     Float_t chisq=0;
@@ -1964,18 +1750,18 @@ void fcnCombiS2(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t if
 
     for (cath=0; cath<2; cath++) {
 //     chisq=0;
-       for (i=0; i<fgNbins[cath]; i++) {
-           Float_t q0=fgCharge[i][cath];
-           Float_t q1=DiscrChargeCombiS2(i,par,cath);
+       for (i=0; i<clusterInput.Nmul(cath); i++) {
+           Float_t q0=clusterInput.Charge(i,cath);
+           Float_t q1=clusterInput.DiscrChargeCombiS2(i,par,cath);
            //      delta=(q0-q1);
            delta=(q0-q1)/q0;
            chisq+=delta*delta;
            qcont+=q1;
            qtot+=q0;
        }
-//     if (cath == 0) chi2temp=chisq/fgNbins[cath];
+//     if (cath == 0) chi2temp=chisq/clusterInput.Nbins[cath];
     }
-//    chisq = chisq/fgNbins[1]+chi2temp;       
+//    chisq = chisq/clusterInput.Nbins[1]+chi2temp;    
     f=chisq;
 }
 
@@ -1985,11 +1771,20 @@ void AliMUONClusterFinderVS::AddRawCluster(const AliMUONRawCluster c)
   // Add a raw cluster copy to the list
   //
     AliMUON *pMUON=(AliMUON*)gAlice->GetModule("MUON");
-    pMUON->AddRawCluster(fChamber,c); 
+    pMUON->AddRawCluster(fInput->Chamber(),c); 
     fNRawClusters++;
     fprintf(stderr,"\nfNRawClusters %d\n",fNRawClusters);
 }
 
+Bool_t AliMUONClusterFinderVS::TestTrack(Int_t t) {
+    if (fTrack[0]==-1 || fTrack[1]==-1) {
+       return kTRUE;
+    } else if (t==fTrack[0] || t==fTrack[1]) {
+       return kTRUE;
+    } else {
+       return kFALSE;
+    }
+}
 
 AliMUONClusterFinderVS& AliMUONClusterFinderVS
 ::operator = (const AliMUONClusterFinderVS& rhs)