loon: AltAlgSliceList.h Source File

00001 
00002 // $Id: AltAlgSliceList.h,v 1.16 2004/11/08 23:22:10 gmieg Exp $
00003 //
00004 // AltAlgSliceList.h
00005 //
00006 //    -- Event Slicing Algorithm
00007 //
00008 // Costas Andreopoulos <C.V.Andreopoulos@rl.ac.uk>
00009 // CCLRC, Rutherford Appleton Laboratory
00010 // July 01, 2003
00012 
00013 #ifndef _ALT_ALG_SLICE_LIST_H_
00014 #define _ALT_ALG_SLICE_LIST_H_
00015 
00016 #include <vector>
00017 #include <map>
00018 
00019 #include "Algorithm/AlgBase.h"
00020 #include "RecoBase/CandStripHandle.h"
00021 #include "MessageService/MsgFormat.h"
00022 
00023 #include <TCanvas.h>
00024 #include <TH1D.h>
00025 
00026 typedef struct centroid { 
00027            double tu, zu, u, tv, zv, v; 
00028 } Centroid_t;
00029 
00030 typedef struct timeslice {  
00031            bool all_sync; 
00032            double tmin, tpeak, tmax, q, qt;  
00033 } TimeSlice_t;
00034 
00035 typedef enum   peakfinderconf { 
00036            kDefault, kLowActivity, kMuSpectrometer 
00037 } PeakFinderConf_t;
00038 
00039 typedef enum wcandslice { 
00040            eNone, eSame, eOther, eMany  
00041 } WhichCandSlice_t;
00042 
00043 typedef enum striptime { 
00044            eTime, eBegTime, eEndTime, eCorrBegTime  
00045 } StripTime_t;
00046 
00047 const int kMuSpec1stPlane = 121;
00048 
00049 class AltAlgSliceList : public AlgBase
00050 {
00051 public:
00052 
00053    AltAlgSliceList();
00054    virtual ~AltAlgSliceList();
00055    virtual void RunAlg(AlgConfig & ac, CandHandle & ch, CandContext & cx);
00056    virtual void Trace(const char * c) const;
00057 //______________________________________________________________________________
00058 
00059 private:
00060 
00061    //-------------------------- methods for recursive time-profile peak finding, 
00062    //                    determination of slice seeds and 0-th order refinement
00063    
00064    std::vector<TimeSlice_t> getSliceSeeds(
00065                               CandStripHandleItr * cshItr, bool recursive_mode, 
00066                                               PeakFinderConf_t conf = kDefault);
00067    void peakFinder(
00068            std::vector<TimeSlice_t> & slice_seeds, PeakFinderConf_t conf) const;
00069    
00070    void reduceTimeWindows(
00071            CandStripHandleItr * cshItr, std::vector<TimeSlice_t> & slice_seeds);
00072 
00073    std::vector<TimeSlice_t> findSmallerPeaks(CandStripHandleItr * cshItr);
00074 
00075    void mergeSeeds(
00076          std::vector<TimeSlice_t> src, std::vector<TimeSlice_t> & target) const;
00077 
00078    void updateSliceSeedInfo(CandStripHandleItr * cshItr, 
00079                                   std::vector<TimeSlice_t> & slice_seeds) const;
00080 
00081    void updateSingleSliceSeedInfo(
00082                          CandStripHandleItr * cshItr, TimeSlice_t & seed) const;
00083 
00084    void removeNullSeeds(std::vector<TimeSlice_t> & slice_seeds) const;
00085 
00086    void resetSliceSeed(TimeSlice_t & seed) const;
00087 
00088    //----------------- method for filling the generic slice-list data structure
00089 
00090    std::map<int, std::vector<CandStripHandle *> > fillSliceSeeds(
00091              CandStripHandleItr cshItr, std::vector<TimeSlice_t> & slice_seeds);
00092 
00093    void sortSlicesInTime(
00094            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00095 
00096    void sortSlicesInPlaneNo(
00097            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00098 
00099    //-------------------------------- methods for time manipulation
00100 
00101    TH1D * createSubsetTimeProfile(CandStripHandleItr * cshItr);
00102 
00103    void getSnarlTimeBoundaries(
00104                CandStripHandleItr * cshItr, double & tmin, double & tmax) const;
00105 
00106    double getStripTime(CandStripHandle * csh, StripTime_t st, 
00107                               StripEnd::StripEnd_t se = StripEnd::kWhole) const;
00108 
00109    //------------------------------- methods for reconstructed slice refinement
00110 
00111    // method for slice refinement via dissolving mini and/or scattered slices
00112    void initSliceFiltering(
00113            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00114 
00115    // method for slice refinement via merging
00116    void sliceMerger(
00117            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00118 
00119    // k-Means - based 3-D clustering for refining strip -> slice assignments
00120    void eventClustering(
00121                  std::map<int, std::vector<CandStripHandle *> > & event_slices);
00122         
00123    // Minimal Spanning Tree based clustering for fine slice-structure search
00124    void sliceSplitter(
00125                  std::map<int, std::vector<CandStripHandle *> > & event_slices);
00126 
00127    //--------------- methods for handling map<int, vector<CandStripHandle *> >
00128 
00129    void Merge(
00130         const std::pair<int, int> & slices_to_merge,
00131            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00132 
00133    void Split(
00134         int slice, std::vector<CandStripHandle *>::size_type n_separator,
00135            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00136 
00137    void Split(
00138         int slice, std::vector< std::vector<CandStripHandle *> > groups,
00139            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00140           
00141    int nextAvailableSliceId(
00142      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00143           
00144    //--------------------------------------- utility methods for slice filtering
00145       
00146    std::vector<CandStripHandle *> dissolveMiniSlices(
00147            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00148 
00149    void giveOrphanStripsForAdoption(
00150                             std::vector<CandStripHandle *> orphan_strips, 
00151            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00152              
00153    int findBestSliceToAdoptAStrip(
00154        bool qweight, int dpln, double dt, CandStripHandle * orphan_strip,
00155            std::map<int, std::vector<CandStripHandle *> > * event_slices) const;
00156 
00157    std::vector<int> * findSlicesToBeEliminated(
00158      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00159 
00160    bool sliceShouldBeEliminated(
00161     std::map<int, std::vector<CandStripHandle *> >::const_iterator slice) const;
00162 
00163    bool smallNumberOfHitStrips(
00164     std::map<int, std::vector<CandStripHandle *> >::const_iterator slice) const;
00165      
00166    bool smallAmountOfCharge(
00167     std::map<int, std::vector<CandStripHandle *> >::const_iterator slice) const;
00168    
00169    bool noStripsInOneView(
00170     std::map<int, std::vector<CandStripHandle *> >::const_iterator slice) const;
00171 
00172    //---------------------------------------- utility methods for slice merging
00173 
00174    std::pair<int, int> findSlicesToMerge(
00175            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00176    
00177    bool checkIfSlicesShouldBeMerged(
00178            std::vector<CandStripHandle *> & slice_0,
00179                                 std::vector<CandStripHandle *> & slice_1) const;
00180         
00181    bool checkSliceDistanceInTimeAndZ(
00182            std::vector<CandStripHandle *> & slice_0,
00183                                 std::vector<CandStripHandle *> & slice_1) const;
00184 
00185    bool checkSliceDistanceInTimeAndSpace(
00186            std::vector<CandStripHandle *> & slice_0,
00187                                 std::vector<CandStripHandle *> & slice_1) const;
00188                                                                 
00189    bool checkPeakTimeDifference(
00190            const std::vector<CandStripHandle *> & slice_0,
00191                           const std::vector<CandStripHandle *> & slice_1) const;
00192    
00193    bool checkPeakZDifference(
00194            const std::vector<CandStripHandle *> & slice_0,
00195                           const std::vector<CandStripHandle *> & slice_1) const;
00196    
00197    bool checkPeakUVDifference(
00198            const std::vector<CandStripHandle *> & slice_0,
00199                           const std::vector<CandStripHandle *> & slice_1) const;
00200 
00201    bool checkLeadingTrailingEdgeDistance(
00202            std::vector<CandStripHandle *> & slice_0,
00203                                 std::vector<CandStripHandle *> & slice_1) const;
00204                                                                 
00205    double averageSliceZ(
00206         std::vector<CandStripHandle *> & slice, double tmin, double tmax) const;
00207 
00208    double minSliceTime(std::vector<CandStripHandle *> & slice) const;
00209    
00210    double maxSliceTime(std::vector<CandStripHandle *> & slice) const;
00211                                 
00212    //---------- utility methods for 3-D *non-hierarchical* (k-Means) clustering 
00213 
00214    bool singleKMeansIteration(
00215                  std::map<int, std::vector<CandStripHandle *> > & event_slices);
00216 
00217    int findBestSliceToHostStrip(
00218                           const pair<int, CandStripHandle *> & slc_strip,
00219                              const std::map<int, Centroid_t> & centroids,
00220      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00221      
00222    std::vector<int> checkForAlternativeSlices(
00223                           const pair<int, CandStripHandle *> & slc_strip,                         
00224      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00225 
00226    WhichCandSlice_t characterizeCandidateSlices(
00227                          int source, const std::vector<int> & alt_slices) const;
00228 
00229    bool needToFilterStrips(
00230      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00231 
00232    bool ZeroChargeInAView(
00233                             const std::vector<CandStripHandle *> & slice) const;     
00234      
00235    bool updateCentroids(std::map<int, Centroid_t> & centroids,
00236                  std::map<int, std::vector<CandStripHandle *> > & event_slices);
00237                                                                   
00238    Centroid_t computeSliceCentroid(
00239                             const std::vector<CandStripHandle *> & slice) const;
00240 
00241    std::map<int, Centroid_t> computeSlicesCentroid(
00242      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00243  
00244    int getClosestCentroid(CandStripHandle * strip,
00245                              const std::map<int, Centroid_t> & centroids) const;
00246       
00247    int selectBestCandidateSlice(
00248                          const pair<int, CandStripHandle *> & slc_strip,  
00249                                    const std::vector<int> & cand_slices,
00250      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00251 
00252    std::multimap<double, int, less<double> > MakeDistanceMap(
00253            CandStripHandle * strip, const std::vector<int> & cand_slices,
00254      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00255            
00256    int moreOccurences(
00257                      std::multimap<double, int> ordered_slice_ids, int n) const;
00258 
00259    void rellocateStrip(int source, int target, CandStripHandle * strip,
00260            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00261 
00262    //------------ utility methods for 3-D *hierarchical* (MST-based) clustering 
00263 
00264    // MST method in my private, non-committed yet,  version. To be 
00265    // committed after the new C++ MC is in place and I can run some checks
00266  
00267    //--------------------------------------- "Near Detector" - specific methods
00268 
00269    void assignMuonSpectrHits2Slices(  CandStripHandleItr cshItr,
00270                  std::map<int, std::vector<CandStripHandle *> > & event_slices);
00271 
00272    std::vector<int> slicesActiveUpstreamOfMuSpec(
00273      const std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00274 
00275    std::map<int, std::vector<CandStripHandle *> > getMuSpecSlices(
00276                                                      CandStripHandleItr cshItr);
00277    void sliceMatcher(
00278              std::map<int, std::vector<CandStripHandle *> > & mu_slices,
00279            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00280 
00281    bool findMatchForMuSpecSlice(
00282     std::vector<CandStripHandle *> mu_strips, std::vector<int> slice_ids, 
00283            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00284 
00285    bool sliceHasStripsInUpstreamMuSpectrometer(
00286                         const std::vector<CandStripHandle *> & mu_strips) const;
00287 
00288    void addMuSpecStripsToExistingSlice( 
00289                         const std::vector<CandStripHandle *> & mu_strips,
00290                            std::vector<CandStripHandle *> & slice_strips) const;
00291 
00292    void newSliceToHostMuSpecStrips( 
00293                         const std::vector<CandStripHandle *> & mu_strips, 
00294            std::map<int, std::vector<CandStripHandle *> > & event_slices) const;
00295 
00296    double averageTimeAtBeginOfDownstreamSlice(
00297                                  std::vector<CandStripHandle *> & strips) const;
00298    
00299    double averageTimeAtEndOfUpstreamSlice(
00300                                  std::vector<CandStripHandle *> & strips) const;
00301 
00302    double averageTime(
00303                       const std::vector<CandStripHandle *> & strips,  
00304                       std::vector<CandStripHandle *>::const_iterator end) const;
00305 
00306    //---------- methods for debugging purposes & sanity tests / mesg formatting
00307 
00308    void printStripList(
00309                    CandStripHandleItr cshItr,  const char * comment = "") const;
00310 
00311    void printSliceSeeds(const std::vector<TimeSlice_t> & slice_seeds, 
00312                         const char * comment = "") const;
00313 
00314    void printSliceSeed(const TimeSlice_t & seed) const;
00315  
00316    void printSlices(
00317      const std::map<int, std::vector<CandStripHandle *> > & event_slices,
00318                                               const char * comments = "") const;
00319    void printSlice(
00320     std::map<int, std::vector<CandStripHandle *> >::const_iterator slice) const;
00321 
00322    void printSlice(const std::vector<CandStripHandle *> & slice) const;
00323 
00324    void printStrip(const CandStripHandle * strip) const;
00325 
00326    void printTimeBin(
00327            int bin, const TH1D * time_profile, const char * comment = "") const;
00328   
00329    std::string asString( WhichCandSlice_t which_slice ) const;
00330    std::string asString( PeakFinderConf_t conf        ) const;
00331 
00332    //------------------ methods for graphics-aided 'online' debugging / testing
00333 
00334    void display(
00335          std::map<int, std::vector<CandStripHandle *> > event_slices,
00336                                                      const char * comment = "");
00337    void eventDisplay(TCanvas * c, 
00338                    std::map<int, std::vector<CandStripHandle *> > event_slices);
00339 
00340    void eventDisplaySingleSlice(TCanvas * c, int sliceKey,
00341                    std::map<int, std::vector<CandStripHandle *> > event_slices);
00342 
00343    void plot3DClusters(TCanvas * c, 
00344                    std::map<int, std::vector<CandStripHandle *> > event_slices);
00345                      
00346    //----------- transform map<int, vector<CandStripHandle *> > into CandSlices
00347 
00348    void buildCandidates(
00349                  std::map<int, std::vector<CandStripHandle *> > & event_slices,
00350                                             CandHandle & ch, CandContext & cx);
00351 
00352    //-------------------------------------------- algorithm config &  self-init
00353 
00354    void getAlgorithmConfiguration(const AlgConfig & ac);   
00355    void internalInit(void);
00356    
00357    //----------------------------------------- algorithm configuration variables   
00358 
00359    int        fGrfxDebugGraphics;
00360    int        fGrfxTimeProfileLogView;
00361    int        fPkfPeakThreshold;
00362    int        fPkfNSuccessiveEmptyBins;
00363    int        fPkfMuSpecPeakThreshold;
00364    int        fPkfMuSpecNSuccessiveEmptyBins;
00365    int        fPkfLowPeakThreshold;
00366    int        fPkfLowNSuccessiveEmptyBins;
00367    int        fPkfNofMergedTimeBins;   
00368    int        fPkfTimeWindowBefPeak;  
00369    int        fPkfTimeWindowAftPeak;  
00370    int        fOrphanStripsPlaneWindow;
00371    int        fKMeansPlaneWindow;
00372    int        fMuSpecNUpstrPlanes;
00373    int        fMuSpecNHitStripsBefSpectr;
00374    int        fMuSpecNPlnBefSpectr;
00375    int        fMinNoHitStrips;
00376    bool       fPkfWeightProfileWithCharge;
00377    bool       fPkfRecursivePeakSearch;
00378    bool       fRefinementDissolving;    
00379    bool       fRefinementMerging;      
00380    bool       fRefinementKMeansClustering;
00381    bool       fRefinementMSTClustering;
00382    bool       fMuSpecSuppressUnmatchedSlices;
00383    bool       fOrphanStripsQWeight;   
00384    double     fOrphanStripsTimeWindow; 
00385    double     fKMeansTimeWindow;
00386    double     fKMeansTPosWindow;      
00387    double     fMuSpecTimeAftUpstrActivity;
00388    double     fMuSpecTimeBefUpstrActivity;
00389    double     fTimeResolution;
00390    double     fTimeDiffBetweenPeaks;
00391    double     fZDiffBetweenPeaks;      
00392    double     fZDiffBetweenEnds;
00393    double     fUVDiffBetweenPeaks;      
00394    double     fMinCharge;         
00395       
00396    //-- internal variables
00397    int        fPeakFinderNestingLevel;  // recursive peak finder nesting level
00398    int        fkMeansIteration;  // k-Means clustering iteration
00399    int        fNTimeBins;
00400    bool       fUpdateCentroids; 
00401    double     fTimeProfileMax;
00402 
00403    TH1D *     fTimeProfile;
00404    double     ftmin;
00405    double     ftmax;   
00406 
00407    TH1D *     fSubsetTimeProfile;
00408    double     fSubsettmin;
00409    double     fSubsettmax;   
00410 
00411    //---------------- message formatting objects   
00412    MsgFormat  fFmtStp;
00413    MsgFormat  fFmtPln;
00414    MsgFormat  fFmtSlc;   
00415    MsgFormat  fFmtTime;
00416    MsgFormat  fFmtQ;
00417 
00418 ClassDef(AltAlgSliceList,0) // AltReco's SliceList Algorithm Class
00419 };
00420 
00421 #endif  // _ALT_ALG_SLICE_LIST_H_