diff --git a/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTKPatternBySector.h b/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTKPatternBySector.h index 64e616a4699839b54e71a9cd5b6225b5db9a7d81..1f98e4ccd3d8dd1b13540968ec3def7ac6ca301d 100644 --- a/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTKPatternBySector.h +++ b/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTKPatternBySector.h @@ -10,16 +10,27 @@ * IO interface for reading and writing patterns to files * For reading patterns, chains of root-files are supported * + * Two basic formats are supported: + * "Forest" there is one TTree per sector -> O(10000) TTrees + * -> excessive buffer memory consumption in root + * "Indexed" there is one TTree with data and two auxillary TTrees + * + * Abstract interface classes: + * FTKPatternBySectorReaderBase basic methods for reading patterns + * FTKPatternBySectorWriterBase basic methods for writing patterns + * use the Factory() methods to create the appropriate + * reader/writer implementation + * + * Implementation: + * FTKPatternBySectorForestReader for reading the "Forest" format + * FTKPatternBySectorIndexedReader for reading the "Indexed" format + * + * FTKPatternBySectorForestWriter for writing the "Forest" format + * FTKPatternBySectorIndexedWriter for writing the "Indexed" format + * * Author: Stefan Schmitt, DESY * Date: 2013/12/02 - * - * classes: - * FTKPatternBySectorBase - * base class for sector-and-coverage ordered patterns root IO - * FTKPatternBySectorReader - * class for reading - * FTKPatternBySectorWriter - * class for writing + * Major modification: 2016/11/18 * */ @@ -54,94 +65,292 @@ protected: CONTENT_TYPE fContentType; }; -class FTKPatternBySectorReader : public FTKPatternBySectorBase { -public: - // typedefs and comparison classes - typedef std::map<int,FTKPatternRootTreeReader *> PatternTreeBySector_t; - class FTKPatternTreeBySectorCompare { - public: - // helper class to compare patterns from different sectors - // and order them by coverage then sector number - explicit FTKPatternTreeBySectorCompare - (FTKHitPatternCompare hitComp=FTKHitPatternCompare()) - : fHitCompare(hitComp) { } - bool operator() - (PatternTreeBySector_t::iterator const &a, - PatternTreeBySector_t::iterator const &b) const { - FTKPatternRootTreeReader const *ra=(*a).second; - FTKPatternRootTreeReader const *rb=(*b).second; - int ca=ra->GetPattern().GetCoverage(); - int cb=rb->GetPattern().GetCoverage(); - // compare coverages - if(ca>cb) return true; - if(ca<cb) return false; - // compare hits - //return fHitCompare(ra->GetPattern().GetHitPattern(), - // rb->GetPattern().GetHitPattern()); //) return true; - // compare sector - return (*a).first<(*b).first; - } - protected: - FTKHitPatternCompare fHitCompare; - }; - typedef std::set<PatternTreeBySector_t::iterator,FTKPatternTreeBySectorCompare> SectorSet_t; - -public: - //explicit FTKPatternBySectorReader(TDirectory &dir, - // char const *name="FTKPatternBySector", - // int nSub=1,int iSub=0); // read pattern bank - explicit FTKPatternBySectorReader(TDirectory &dir, - char const *name="FTKPatternBySector", - std::set<int> const *sectorList=0); // read pattern bank - explicit FTKPatternBySectorReader(FTKRootFileChain &chain); // read pattern bank - ~FTKPatternBySectorReader(); - - void SelectSubregion(int nSub,int iSub); - - Long64_t GetNPatterns(void) const; - Long64_t GetNPatterns(int sector) const; - void GetNPatternsByCoverage(std::map<int,int> &coverageMap) const; - void GetNPatternsBySectorCoverage(std::map<int,std::map<int,int> >§orCoverageMap) const; +class FTKPatternBySectorReader : public virtual FTKPatternBySectorBase { + public: + static FTKPatternBySectorReader *Factory + (TDirectory &dir,std::set<int> const *sectorList=0); + static FTKPatternBySectorReader *Factory(FTKRootFileChain &chain); + + virtual ~FTKPatternBySectorReader(); + virtual Long64_t GetNPatternsTotal(void) const; + + // read patterns (optimized IO), given a minimum number of patterns + // repeated calls to this function will result in all sectors to be read + // for sectors appended to the list, all patterns have been read + // (use Rewind() to start all over) + // returns: next sector to read (or -1) + virtual int ReadRaw(int firstSector, + std::list<FTKPatternOneSector *> §orList, + uint64_t maxPattern); + + // for the given sector, read all patterns down to minCoverage + // minCoverage<=1: read all data (also works for non-merged banks) + virtual FTKPatternOneSector *Read(int sector,int minCoverage)=0; + + // read all patterns of a sector, also works for non-merged banks + FTKPatternOneSector *MergePatterns(int sector) { + RewindSector(sector); return Read(sector,0); + } + + // rewind one/all sectors + void Rewind(void); + virtual void RewindSector(int sector)=0; + // determine multiplicity of each coverage + void GetNPatternsByCoverage(std::map<int,int> &coverageMap) const; + // for each coverage, list sectors and number of patterns + typedef std::map<int,int> SectorMultiplicity_t; + typedef std::map<int,SectorMultiplicity_t> SectorByCoverage_t; + virtual SectorByCoverage_t const &GetSectorByCoverage(void) const { + return m_SectorByCoverage; + } + // export patterns as ASCII int WritePatternsToASCIIstream(std::ostream &out, int iSub,int nSub); // export as ASCII - // save data from one sector, ordered by coverage - FTKPatternOneSector *MergePatterns(int sector) const; // get data from one sector - // the following methods are for looping over the data of one sector - int GetFirstSector(void) const; - int GetNextSector(int sector) const; - Long64_t GetPatternByCoverage(SectorSet_t& sectorByCoverage, int iSub, int nSub) ; + // total number of patterns in a sector + virtual Long64_t GetNPatterns(int sector) const=0; + // for loop over all sectors + virtual int GetFirstSector(void) const=0; + virtual int GetNextSector(int sector) const=0; -protected: - PatternTreeBySector_t fPatterns; + // for loops over coverage, within coverage over sector + int GetMaxCoverage(void) const; + int GetNextCoverage(int coverage) const; + int GetCoverageFirstSector(int coverage) const; + int GetCoverageNextSector(int coverage,int sector) const; + protected: + FTKPatternBySectorReader(char const *name); + SectorByCoverage_t m_SectorByCoverage; +}; +class FTKPatternBySectorWriter : public virtual FTKPatternBySectorBase { + public: + enum WRITERTYPE_t { + kWriterDefault=0, + kWriterForest=1, + kWriterIndexed=2 + }; + static FTKPatternBySectorWriter *Factory + (TDirectory &dir,WRITERTYPE_t type=kWriterDefault); + // append patterns from file (down to minCoverage) + // the default algorithm merges one sector a time, calling + // AppendMergedPatterns() fro each sector. + // This procedure is not optimal for all file formats + // For this reason the method is virtual, to allow for + // a more efficient implementation + virtual int AppendMergedPatterns + (FTKPatternBySectorReader &in,int minCoverage); + // append patterns from ASCII file + int AppendPatternsFromASCIIstream(std::istream &in); + // append a single pattern + int AppendPattern(FTKPatternWithSector const &pattern) { + return AppendPattern(pattern.GetSectorID(),pattern);} + int AppendPattern(int sector,FTKPatternWithCoverage const &pattern) { + return AppendPattern(sector,pattern.GetCoverage(), + pattern.GetHitPattern()); + } + virtual int AppendPattern(int sector,int coverage, + FTKHitPattern const &pattern)=0; + protected: + FTKPatternBySectorWriter(const char *name,TDirectory &dir) + : FTKPatternBySectorBase(name),fDir(dir) { } + virtual int AppendMergedPatternsSector + (int sector,FTKPatternOneSectorOrdered const *ordered); + TDirectory &fDir; + static uint64_t PATTERN_CHUNK; }; -class FTKPatternBySectorBlockReader : public FTKPatternBySectorReader { +class FTKPatternBySectorForestWriter : public FTKPatternBySectorWriter { +public: + explicit FTKPatternBySectorForestWriter(TDirectory &dir); // create pattern bank + virtual ~FTKPatternBySectorForestWriter(); // save patterns to TDirectory + virtual int AppendPattern(int sector,int coverage, + FTKHitPattern const &pattern); + protected: + virtual int AppendMergedPatternsSector(int sector, + FTKPatternOneSectorOrdered const *ordered); + typedef std::map<int,FTKPatternRootTree *> PatternTreeBySectorRW_t; + PatternTreeBySectorRW_t fPatterns; +}; + +class FTKPatternBySectorForestReader : public FTKPatternBySectorReader { public: - FTKPatternBySectorBlockReader - (TDirectory &dir,std::set<int> const *sectorList); - void Rewind(void); - FTKPatternOneSector *Read(int sector,int minCoverage); + // for this reader, the sectors are organized in a forest of TTrees + explicit FTKPatternBySectorForestReader + (TDirectory &dir,std::set<int> const *sectorList=0,int *error=0); // read pattern bank + explicit FTKPatternBySectorForestReader(FTKRootFileChain &chain, + int *error=0); // read pattern bank + virtual ~FTKPatternBySectorForestReader(); + virtual Long64_t GetNPatterns(int sector) const; + virtual int GetFirstSector(void) const; + virtual int GetNextSector(int sector) const; + virtual void RewindSector(int sector=-1); + virtual FTKPatternOneSector *Read(int sector,int minCoverage); + protected: + typedef std::map<int,FTKPatternRootTreeReader *> PatternTreeBySectorRO_t; + PatternTreeBySectorRO_t fPatterns; + void InitializeSectorByCoverageTable(void); }; -class FTKPatternBySectorWriter : public FTKPatternBySectorBase { +class FTKPatternBySectorIndexed : public virtual FTKPatternBySectorBase { + public: + //int GetBitsPerPlane(void) const; + static char const *s_cdataTreeName; + static char const *s_indexTreeName; + static char const *s_cnpattBranchName; + static char const *s_cnpattBranchID; + static char const *GetCPatternBranchName(uint32_t plane); + static char const *GetCPatternBranchID(uint32_t plane); + static char const *s_sectorBranchName; + static char const *s_sectorBranchID; + static char const *s_ncovBranchName; + static char const *s_ncovBranchID; + static char const *s_npatternBranchName; + static char const *s_npatternBranchID; + static char const *s_coverageBranchName; + static char const *s_coverageBranchID; + static char const *GetNDecoderBranchName(uint32_t plane); + static char const *GetNDecoderBranchID(uint32_t plane); + static char const *GetDecoderDataBranchName(uint32_t plane); + static char const *GetDecoderDataBranchID(uint32_t plane); + protected: + static std::vector<std::string> s_cpatternBranch,s_cpatternID; + static std::vector<std::string> s_ndecoderBranch,s_ndecoderID; + static std::vector<std::string> s_decoderDataBranch,s_decoderDataID; +}; + +class FTKPatternBySectorIndexedWriter : public FTKPatternBySectorWriter, + FTKPatternBySectorIndexed { public: - explicit FTKPatternBySectorWriter(TDirectory &dir); // read or create pattern bank - virtual ~FTKPatternBySectorWriter(); // save patterns to TDirectory - int AppendMergedPatterns(FTKPatternBySectorReader &in,int minCoverage); // append patterns - int AppendPatternsFromASCIIstream(std::istream &in); // append - int AppendMergedPatterns(int sector, - FTKPatternOneSectorOrdered const *ordered); - int AppendPattern(int sector,FTKPatternWithCoverage const &pattern); - int AppendPattern(FTKPatternWithSector const &pattern) { return AppendPattern(pattern.GetSectorID(),pattern);} - TDirectory& GetTDirectoryAddress() { return fDir;} - -protected: - TDirectory &fDir; - typedef std::map<int,FTKPatternRootTree *> PatternTreeBySector_t; - PatternTreeBySector_t fPatterns; + explicit FTKPatternBySectorIndexedWriter(TDirectory &dir); // create pattern bank + virtual ~FTKPatternBySectorIndexedWriter(); // save patterns to TDirectory + virtual int AppendPattern(int sector,int coverage, + FTKHitPattern const &pattern); + protected: + TTree *m_cpatternDataTree; + std::vector<TBranch *> m_cpatternBranch; + TBranch *m_cnpattBranch; + + TTree *m_cpatternIndexTree; + TBranch *m_npatternBranch; + TBranch *m_coverageBranch; + std::vector<TBranch*> m_decoderDataBranch; + // total number of patterns in memory + class CodedPatternStorage { + public: + explicit CodedPatternStorage(uint32_t nPlane) : m_data(nPlane) { + for(size_t plane=0;plane<m_data.size();plane++) + m_data[plane].resize(15000000); + m_next=0; + } + // returns size of storage + uint32_t GetSize(void) const { + return m_data[0].size(); + } + + // returns 0 is storage is full + // index in stoage otherwise + inline uint32_t AllocatePattern(void) { + uint32_t r=(m_next+1)%(m_data[0].size()+1); + if(r) m_next=r; + return r; + } + inline void Store(uint32_t ipatt,uint32_t iplane,uint8_t data) { + m_data[iplane][ipatt-1]=data; + } + inline uint8_t Get(uint32_t ipatt,uint32_t iplane) const { + return m_data[iplane][ipatt-1]; + } + inline void Delete(uint32_t ipatt) { + if(ipatt && (ipatt==m_next)) --m_next; + } + inline void Clear(void) { m_next=0; } + bool operator()(uint32_t a,uint32_t b) const { + for(size_t plane=0;plane<m_data.size();plane++) { + if(Get(a,plane)<Get(b,plane)) return true; + if(Get(a,plane)>Get(b,plane)) return false; + } + return false; + } + protected: + uint32_t m_next; + std::vector<std::vector<uint8_t> > m_data; + }; + CodedPatternStorage *m_patternStorage; + typedef std::map<uint32_t,uint32_t,CodedPatternStorage&> + CPatternWithCoverage_t; + // indexing + typedef struct PatternData { + PatternData(int nLayer,CodedPatternStorage *storage) + : m_encoder(nLayer),m_data(*storage) { + for(size_t i=0;i<m_encoder.size();i++) { + m_encoder[i].first.resize(1<<8); + } + } + // encoder, translates [SSID] -> uint8_t + std::vector< std::pair<std::vector<int>,std::map<int,int> > > m_encoder; + // patterns with coverage + CPatternWithCoverage_t m_data; + } PatternData_t; + // PatternData indexed by sector + std::map<int,PatternData_t> m_patternData; + uint32_t m_treeSector; + uint32_t m_treeNCov; + std::vector<uint32_t> m_treeCoverage; + std::vector<uint32_t> m_treeNPattern; + std::vector<uint32_t> m_decoderSize; + void Flush(void); +}; + +class FTKPatternBySectorIndexedReader : public FTKPatternBySectorReader, + FTKPatternBySectorIndexed { + public: + // for this reader, the sectors are organized in two TTrees + explicit FTKPatternBySectorIndexedReader + (TDirectory &dir,std::set<int> const *sectorList=0,int *error=0); // read pattern bank + explicit FTKPatternBySectorIndexedReader(FTKRootFileChain &chain, + int *error=0); // read pattern bank + virtual ~FTKPatternBySectorIndexedReader(); + virtual Long64_t GetNPatterns(int sector) const; + virtual int GetFirstSector(void) const; + virtual int GetNextSector(int sector) const; + virtual void RewindSector(int sector=-1); + virtual FTKPatternOneSector *Read(int sector,int minCoverage); + // optimized reader of multiple sectors + virtual int ReadRaw(int firstSector, + std::list<FTKPatternOneSector *> §orList, + uint64_t maxPattern); + protected: + TTree *m_dataTree; + int InitializeIndex(TTree *indexTree,std::set<int> const *sectorList); + + // + std::vector<uint32_t> m_decoderData; + std::vector<uint32_t *> m_decoderDataPtr; + struct DataBlock_t { + // pointer to decoder table + // decoded_SSID = m_decoderDataPtr[plane][coded_SSID] + uint32_t **m_decoderDataPtr; + // pointer to offset in data TTree + uint32_t m_cpatternEntry; + // number of patterns in this Block + uint32_t m_nPattern; + }; + // first index: sector + // second index: coverage (multiple entries are possible - decoder change) + + typedef std::multimap<uint32_t,DataBlock_t> CoverageTable_t; + typedef struct SectorInfo { + // table ordered by coverage + CoverageTable_t m_coverageTable; + // pointer to one entry of the table above (for reading) + CoverageTable_t::const_reverse_iterator m_coveragePtr; + uint64_t m_nPattern; + } SectorInfo_t; + typedef std::map<uint32_t,SectorInfo_t> SectorTable_t; + SectorTable_t m_sectorTable; + FTKPatternWithCoverage *m_patternData; }; #endif diff --git a/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTKPatternOneSector.h b/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTKPatternOneSector.h index 181c31632c1d490e4c64436605bce202f1f45068..5a3655d702f1c0ee10adf9314ab404c2333fc0ac 100644 --- a/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTKPatternOneSector.h +++ b/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTKPatternOneSector.h @@ -65,7 +65,7 @@ public: return (**p).first; } inline int GetCoverage(Ptr_t const &p) const { return (**p).second; } - + inline int GetNumberOfPatterns(void) const { return fPatterns.size(); } int GetSummedCoverage(void) const; protected: OrderedPatternSet_t fPatterns; diff --git a/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTK_CompressedAMBank.h b/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTK_CompressedAMBank.h index 510906f847900d070bf11f3c9de473a5c8c6bfcf..fc26a2cd36d905beb2ec5324748ebad2dc0ce5da 100644 --- a/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTK_CompressedAMBank.h +++ b/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/FTK_CompressedAMBank.h @@ -31,6 +31,7 @@ class TSPMap; class FTKHitPattern; class FTKHitPatternCompare; class FTKPatternOneSector; +class FTKPatternBySectorReader; //#include <boost/container/map.hpp> //#include <boost/container/vector.hpp> diff --git a/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/tsp/TSPROOTBankGenerator.h b/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/tsp/TSPROOTBankGenerator.h index 987e2dce95c33ca6bc3b9f039c16d34307768a12..160796ad50cb14cd9b36f270dd22c083e5714b04 100644 --- a/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/tsp/TSPROOTBankGenerator.h +++ b/Trigger/TrigFTK/TrigFTKSim/TrigFTKSim/tsp/TSPROOTBankGenerator.h @@ -78,7 +78,7 @@ private: std::stringstream m_patternStringID; PatternMap m_patternMap; std::vector<int> m_ssPattern; - FTKPatternBySectorReader* m_preader; + FTKPatternBySectorForestReader* m_preader; }; } diff --git a/Trigger/TrigFTK/TrigFTKSim/scripts/TrigFTKSim_tf.py b/Trigger/TrigFTK/TrigFTKSim/scripts/TrigFTKSim_tf.py index 8c4b1ef98930121405302b95ac24d347ed8ed762..de82fcfe1cf434c7100257ca31a7bb4a00f610cd 100755 --- a/Trigger/TrigFTK/TrigFTKSim/scripts/TrigFTKSim_tf.py +++ b/Trigger/TrigFTK/TrigFTKSim/scripts/TrigFTKSim_tf.py @@ -3,7 +3,7 @@ # Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration ## FTK Simulation Transform -# @version $Id: TrigFTKSim_tf.py 770308 2016-08-26 12:32:36Z vcavalie $ +# @version $Id: TrigFTKSim_tf.py 788284 2016-12-06 19:04:09Z sschmitt $ import argparse import sys @@ -102,6 +102,8 @@ def addFTKSimulationArgs(parser): help='CachePath', group='TrigFTKSim') parser.add_argument("--CachedBank", type=trfArgClasses.argFactory(trfArgClasses.argBool, runarg=True), help="Interpret the pattern bank has a cache", group="TrigFTKSim") + parser.add_argument("--FixEndCapL0", type=trfArgClasses.argFactory(trfArgClasses.argBool, runarg=True), + help="Fix problem with pixel layer assignment in endcap", group="TrigFTKSim") parser.add_argument("--SectorsAsPatterns", type=trfArgClasses.argFactory(trfArgClasses.argInt, runarg=True), help="If 1 allows to use a list of sectors as pattern bank, default 0") diff --git a/Trigger/TrigFTK/TrigFTKSim/share/skeleton.BS_FTK_Creator.py b/Trigger/TrigFTK/TrigFTKSim/share/skeleton.BS_FTK_Creator.py index de5af3da08c3f03edb152d5592259a23c2a0ff22..252856bcd68f748b20572c9d779c044ca69f1396 100755 --- a/Trigger/TrigFTK/TrigFTKSim/share/skeleton.BS_FTK_Creator.py +++ b/Trigger/TrigFTK/TrigFTKSim/share/skeleton.BS_FTK_Creator.py @@ -22,6 +22,18 @@ else: ftkLog.info("Running on all the events") athenaCommonFlags.EvtMax = -1 +## Pre-exec +if hasattr(runArgs,"preExec"): + ftkLog.info("transform pre-exec") + for cmd in runArgs.preExec: + ftkLog.info(cmd) + exec(cmd) + +## Pre-include +if hasattr(runArgs,"preInclude"): + for fragment in runArgs.preInclude: + include(fragment) + if hasattr(runArgs, "skipEvents"): athenaCommonFlags.SkipEvents.set_Value_and_Lock(runArgs.skipEvents) elif hasattr(runArgs, "firstEvent"): diff --git a/Trigger/TrigFTK/TrigFTKSim/share/skeleton.RDO_FTK_Creator.py b/Trigger/TrigFTK/TrigFTKSim/share/skeleton.RDO_FTK_Creator.py index f4c71321af270e0cbfd6451f6c97836dd1e15df3..6c0fa20a4c11573a0a1e5d95d0bff276c5613b65 100755 --- a/Trigger/TrigFTK/TrigFTKSim/share/skeleton.RDO_FTK_Creator.py +++ b/Trigger/TrigFTK/TrigFTKSim/share/skeleton.RDO_FTK_Creator.py @@ -9,7 +9,13 @@ pmjp.PerfMonFlags.doSemiDetailedMonitoring = True from AthenaCommon.Logging import logging ftkLog = logging.getLogger('FTKRDOCreator') #ftkLog.propagate = False -ftkLog.info( '********** STARTING FTKStandaloneMerge **********' ) +ftkLog.info( '********** STARTING FTKStandaloneMerge RDO_FTK **********' ) + +from AthenaCommon.DetFlags import DetFlags +DetFlags.makeRIO.pixel_setOn() +DetFlags.makeRIO.SCT_setOn() +DetFlags.detdescr.all_setOn() +DetFlags.geometry.all_setOn() athenaCommonFlags.FilesInput = runArgs.inputRDOFile @@ -20,6 +26,18 @@ else: ftkLog.info("Running on all the events") athenaCommonFlags.EvtMax = -1 +## Pre-exec +if hasattr(runArgs,"preExec"): + ftkLog.info("transform pre-exec") + for cmd in runArgs.preExec: + ftkLog.info(cmd) + exec(cmd) + +## Pre-include +if hasattr(runArgs,"preInclude"): + for fragment in runArgs.preInclude: + include(fragment) + if hasattr(runArgs, "skipEvents"): athenaCommonFlags.SkipEvents.set_Value_and_Lock(runArgs.skipEvents) elif hasattr(runArgs, "firstEvent"): diff --git a/Trigger/TrigFTK/TrigFTKSim/src/FTKMergeRoot.cxx b/Trigger/TrigFTK/TrigFTKSim/src/FTKMergeRoot.cxx index b941f9ad8785a86838397385020fbacadd601027..c6c144757cd3dc999c75569fbc107e8152bf975d 100644 --- a/Trigger/TrigFTK/TrigFTKSim/src/FTKMergeRoot.cxx +++ b/Trigger/TrigFTK/TrigFTKSim/src/FTKMergeRoot.cxx @@ -128,7 +128,7 @@ int FTKMergeRoot::DoMerge(int MinCoverage,int compression){ } // Sector Reader - FTKPatternBySectorReader* input = new FTKPatternBySectorReader(chain); + FTKPatternBySectorReader* input = FTKPatternBySectorReader::Factory(chain); if(!input->GetNLayers()) { // sanity check Warning("DoMerge")<<"number of layers not set, cannot merge patterns!\n"; delete input; @@ -176,7 +176,7 @@ int FTKMergeRoot::DoMerge(int MinCoverage,int compression){ } Info("DoMerge")<<"merge patterns to root file "<<outFileName<<"\n"; FTKPatternBySectorWriter* patternWriter = - new FTKPatternBySectorWriter(*rootOutputFile); + FTKPatternBySectorWriter::Factory(*rootOutputFile); patternWriter->AppendMergedPatterns(*input,MinCoverage); delete patternWriter; // write to disk patternWriter=0; @@ -270,7 +270,7 @@ void FTKMergeRoot::DoTextImport(){ <<"\" for writing imported patterns\n"; return; } - patternWriter = new FTKPatternBySectorWriter(*tmpFile); + patternWriter = FTKPatternBySectorWriter::Factory(*tmpFile); } Info("DoTextImport") @@ -337,7 +337,7 @@ void FTKMergeRoot::DoTextExport(std::string const &TextOutFilename, int /*MinCov chain.AddFile(inputFileName); } } - FTKPatternBySectorReader* input = new FTKPatternBySectorReader(chain); + FTKPatternBySectorReader* input = FTKPatternBySectorReader::Factory(chain); if(!input) { Fatal("DoTextExport")<<"do not know where to read the pattern from\n"; return; diff --git a/Trigger/TrigFTK/TrigFTKSim/src/FTKPatternBySector.cxx b/Trigger/TrigFTK/TrigFTKSim/src/FTKPatternBySector.cxx index 8d2d31714382c3231c542767e7d94f797b6a90f1..5114fcf660b3b96263ab842ac4b0f8580ac5a746 100644 --- a/Trigger/TrigFTK/TrigFTKSim/src/FTKPatternBySector.cxx +++ b/Trigger/TrigFTK/TrigFTKSim/src/FTKPatternBySector.cxx @@ -16,32 +16,13 @@ #include "TrigFTKSim/FTKRootFile.h" #include <TList.h> #include <TFile.h> +#include <TChain.h> #include <TObjString.h> - -#ifdef DIAG_VMEM #include <sys/time.h> #include <sys/resource.h> -static void printVmemUsage() { - struct rusage usage; - getrusage(RUSAGE_SELF,&usage); - std::cout<<"VMEM usage:" - <<" "<<usage.ru_maxrss/1024./1024.<<"G" - <<"\n"; -} -#endif using namespace std; -static TDirectory *OpenRootFileReadonly(char const *path) { - TFile *r=new TFile(path,"READ"); - if(r && r->IsZombie()) { - // if TFile can not be opened, return 0 - delete r; - r=0; - } - return r; -} - //================== class FTKPatternBySectorBase ======================== FTKPatternBySectorBase::~FTKPatternBySectorBase() { @@ -76,11 +57,433 @@ int FTKPatternBySectorBase::SetNLayers(int nLayer) { return error; } -//================== class FTKPatternBySectorReader ======================== +//================== class FTKPatternBySectorReader ==== + +FTKPatternBySectorReader::FTKPatternBySectorReader(char const *name) + : FTKPatternBySectorBase(name) { +} + +FTKPatternBySectorReader::~FTKPatternBySectorReader() { +} + +Long64_t FTKPatternBySectorReader::GetNPatternsTotal(void) const { + // return total number of patterns + Long64_t n=0; + for(int sector=GetFirstSector();sector>=0;sector=GetNextSector(sector)) { + n += GetNPatterns(sector); + } + return n; +} + +void FTKPatternBySectorReader::GetNPatternsByCoverage +(map<int,int> &coverageMap) const { + coverageMap.clear(); + for(SectorByCoverage_t::const_iterator c=m_SectorByCoverage.begin(); + c!=m_SectorByCoverage.end();c++) { + int &n=coverageMap[(*c).first]; + for(SectorMultiplicity_t::const_iterator s=(*c).second.begin(); + s!=(*c).second.end();s++) { + n+= (*s).second; + } + } +} + +int FTKPatternBySectorReader::GetMaxCoverage(void) const { + int r=0; + if(!m_SectorByCoverage.empty()) + r=(*m_SectorByCoverage.rbegin()).first; + return r; +} + +int FTKPatternBySectorReader::GetNextCoverage(int coverage) const { + SectorByCoverage_t::const_iterator c=m_SectorByCoverage.find(coverage); + int r=0; + if(c!=m_SectorByCoverage.end()) { + --c; + if(c!=m_SectorByCoverage.end()) r=(*c).first; + } + return r; +} + +int FTKPatternBySectorReader::GetCoverageFirstSector(int coverage) const { + int r=-1; + SectorByCoverage_t::const_iterator c=m_SectorByCoverage.find(coverage); + if(c!=m_SectorByCoverage.end()) { + SectorMultiplicity_t::const_iterator s=(*c).second.begin(); + if(s!=(*c).second.end()) r=(*s).first; + } + return r; +} + +int FTKPatternBySectorReader::GetCoverageNextSector +(int coverage,int sector) const { + int r=-1; + SectorByCoverage_t::const_iterator c=m_SectorByCoverage.find(coverage); + if(c!=m_SectorByCoverage.end()) { + SectorMultiplicity_t::const_iterator s=(*c).second.find(sector); + if(s!=(*c).second.end()) { + ++s; + if(s!=(*c).second.end()) { + r=(*s).first; + } + } + } + return r; +} + +int FTKPatternBySectorReader::WritePatternsToASCIIstream +(ostream &out,int iSub,int nSub) { + // read all patterns from all sectors + // write data ordered by + // coverage + // then sector number + // then keep the order of the hits in the input TTree + int error=0; + if(GetContentType()==CONTENT_NOTMERGED) { + Error("WritePatternsToASCIIstream") + <<"data can not be saved as ASCII (still have to be merged)\n"; + } + // write header + out<<GetNPatternsTotal()<<" "<<GetNLayers()<<"\n"; + + // write patterns + int number=0; + Rewind(); + // coverage loop + for(int coverage=GetMaxCoverage();coverage; + coverage=GetNextCoverage(coverage)) { + // sector loop + for(int sector=GetCoverageFirstSector(coverage);sector>=0; + sector=GetCoverageNextSector(coverage,sector)) { + if((nSub>1)&&(sector%nSub)!=iSub) continue; + FTKPatternOneSector *data=Read(sector,coverage); + // pattern loop + for(FTKPatternOneSector::Ptr_t i=data->Begin();i!=data->End();i++) { + // progress monitor + if(!(number %1000)) { + ShowProgress(TString::Format("%d/%d/%d",number,sector,coverage)); + } + // write pattern + out<<number; + FTKHitPattern const &patternData=data->GetHitPattern(i); + for(unsigned iLayer=0;iLayer<patternData.GetSize();iLayer++) { + out<<' '<<patternData.GetHit(iLayer); + } + out<<' '<<sector<<' '<<coverage<<"\n"; + number++; + } + delete data; + if(out.fail()) { + error++; + Error("WritePatternsToASCIIstream")<<"write operation failed\n"; + break; + } + } + } + Info("WritePatternsToASCIIstream") + <<" total patterns written: "<<number<<"\n"; + return error; +} + +void FTKPatternBySectorReader::Rewind(void) { + for(int sector=GetFirstSector();sector>=0;sector=GetNextSector(sector)) { + RewindSector(sector); + } +} + +FTKPatternBySectorReader *FTKPatternBySectorReader::Factory +(TDirectory &dir,set<int> const *sectorList) { + FTKPatternBySectorReader *reader; + int errorIndex=0; + int errorForest=0; + reader=new FTKPatternBySectorIndexedReader(dir,sectorList,&errorIndex); + if(errorIndex) { + delete reader; + reader=new FTKPatternBySectorForestReader(dir,sectorList,&errorForest); + if(errorForest) { + reader=new FTKPatternBySectorIndexedReader(dir,sectorList,0); + delete reader; + reader=new FTKPatternBySectorForestReader(dir,sectorList,&errorForest); + delete reader; + reader=0; + } + } + return reader; +} + +FTKPatternBySectorReader *FTKPatternBySectorReader::Factory +(FTKRootFileChain &chain) { + FTKPatternBySectorReader *reader; + int errorIndex=0; + int errorForest=0; + reader=new FTKPatternBySectorIndexedReader(chain,&errorIndex); + if(errorIndex) { + delete reader; + reader=new FTKPatternBySectorForestReader(chain,&errorForest); + if(errorForest) { + ostream &out=reader->Error("FTKPatternBySectorReader"); + out<<"Can not read some files of the chain ["; + for(int i=0;i<chain.GetLength();i++) { + if(i) out<<","; + out<<chain[i]; + } + out<<"]\n"; + reader=new FTKPatternBySectorIndexedReader(chain,0); + delete reader; + reader=new FTKPatternBySectorForestReader(chain,0); + delete reader; + reader=0; + } + } + return reader; +} + +int FTKPatternBySectorReader::ReadRaw +(int firstSector,list<FTKPatternOneSector *> §orList,uint64_t minPattern) { + int sector; + uint64_t total=0; + for(sector=firstSector;sector>=0;sector=GetNextSector(sector)) { + FTKPatternOneSector *sectorData=Read(sector,0); + if(!sectorData) { + sector=-1; + break; + } + sectorList.push_back(sectorData); + total+= sectorData->GetNumberOfPatterns(); + if(total>=minPattern) break; + } + return sector; +} + + +//================== class FTKPatternBySectorWriter ==== + +FTKPatternBySectorWriter *FTKPatternBySectorWriter::Factory +(TDirectory &dir,WRITERTYPE_t type) { + if(type==kWriterForest) + return new FTKPatternBySectorForestWriter(dir); + else if(type==kWriterIndexed) + return new FTKPatternBySectorIndexedWriter(dir); + // default + return new FTKPatternBySectorIndexedWriter(dir); +} + +#ifdef OLD_ALGO +int FTKPatternBySectorWriter::AppendMergedPatterns +(FTKPatternBySectorReader &source,int minCoverage) { + // loop over all sectors + // for each sector, import-and-merge patterns + int error=0; + SetNLayers(source.GetNLayers()); + //Info("AppendMergedPatterns")<<"loop over sectors\n"; + int nTotal=0; + int64_t totalCoverage=0; + int printed=0; + double nextPrint=1.0; + static int printSector=0; + + struct rusage startUsage; + int numSectorUsage=0; + getrusage(RUSAGE_SELF,&startUsage); + + for(int sector=source.GetFirstSector();sector>=0; + sector=source.GetNextSector(sector)) { + //ShowProgress(TString::Format("sector %d",sector)); + FTKPatternOneSector *mergedPatterns=source.MergePatterns(sector); + if(printSector) { + printSector--; + int n=mergedPatterns->GetNumberOfPatterns(); + FTKPatternOneSector::Ptr_t p=mergedPatterns->Begin(); + for(int k=0;k<2;k++) { + for(int i=0;(i<5)&&(p!=mergedPatterns->End());i++,p++) { + cout<<mergedPatterns->GetNumberOfPatterns()-n + <<" "<<sector<<" "<<mergedPatterns->GetCoverage(p); + for(int j=0;j<GetNLayers();j++) { + cout<<" "<<mergedPatterns->GetHitPattern(p).GetHit(j); + } + cout<<"\n"; + } + --n; + if(n>10) { + cout<<" ...\n"; + while(n>10) { + --n; + ++p; + } + } + } + cout<<" ...\n"; + } + nTotal += mergedPatterns->GetNumberOfPatterns(); + FTKPatternOneSectorOrdered *orderedByCoverage= + mergedPatterns->OrderPatterns(FTKPatternOrderByCoverage(minCoverage)); + totalCoverage += orderedByCoverage->GetSummedCoverage(); + printed++; + struct rusage currentUsage; + numSectorUsage++; + getrusage(RUSAGE_SELF,¤tUsage); + struct timeval tdiff; + timersub(¤tUsage.ru_utime,&startUsage.ru_utime,&tdiff); + double dt=(tdiff.tv_sec+tdiff.tv_usec/1.E6)/numSectorUsage; + if(printed>=nextPrint) { + ShowProgress(TString::Format + ("sector %d patterns %d total %d %ld time/sector %.3f (nsector=%d)", + sector, + mergedPatterns->GetNumberOfPatterns(),nTotal, + totalCoverage,dt,numSectorUsage)); + nextPrint=printed*1.3; + } + error=AppendMergedPatternsSector(sector,orderedByCoverage); + if(error) {delete orderedByCoverage; delete mergedPatterns; break;} + delete orderedByCoverage; + delete mergedPatterns; + } + Info("AppendMergedPatterns") + <<"nsector="<<printed<<" npattern="<<nTotal + <<" npattern*coverage="<<totalCoverage<<"\n"; + // save sector index as TTree + return error; +} +#else + +uint64_t FTKPatternBySectorWriter::PATTERN_CHUNK=20000000; + + +int FTKPatternBySectorWriter::AppendMergedPatterns +(FTKPatternBySectorReader &source,int minCoverage) { + // loop over all sectors + // for each sector, import-and-merge patterns + int error=0; + SetNLayers(source.GetNLayers()); + Info("AppendMergedPatterns")<<"loop over sectors\n"; + int nTotal=0; + int64_t totalCoverage=0; + int printed=1; + double nextPrint=1.0; + struct rusage startUsage; + int numSectorUsage=0; + getrusage(RUSAGE_SELF,&startUsage); + + int readSector=source.GetFirstSector(); + while(readSector>=0) { + std::list<FTKPatternOneSector *> sectorList; + readSector=source.ReadRaw(readSector,sectorList,PATTERN_CHUNK); + int nSector=0; + int nPattern=0; + for(std::list<FTKPatternOneSector *>::iterator iSector=sectorList.begin(); + iSector!=sectorList.end();iSector++) { + FTKPatternOneSector *mergedPatterns=(*iSector); + nTotal += mergedPatterns->GetNumberOfPatterns(); + FTKPatternOneSectorOrdered *orderedByCoverage= + mergedPatterns->OrderPatterns + (FTKPatternOrderByCoverage(minCoverage)); + totalCoverage += orderedByCoverage->GetSummedCoverage(); + error+=AppendMergedPatternsSector(mergedPatterns->GetSector(), + orderedByCoverage); + delete orderedByCoverage; + delete (*iSector); + nSector++; + nPattern+= mergedPatterns->GetNumberOfPatterns(); + } + struct rusage currentUsage; + numSectorUsage+=nSector; + getrusage(RUSAGE_SELF,¤tUsage); + struct timeval tdiff; + timersub(¤tUsage.ru_utime,&startUsage.ru_utime,&tdiff); + double dt=(tdiff.tv_sec+tdiff.tv_usec/1.E6)/numSectorUsage; + if(printed>=nextPrint) { + ShowProgress(TString::Format + ("patterns %d total %d %ld time/sector %.3f (nsector=%d)", + nPattern,nTotal, + totalCoverage,dt,numSectorUsage)); + nextPrint=printed*1.3; + } + printed++; + } + Info("AppendMergedPatterns") + <<"nsector="<<printed<<" npattern="<<nTotal + <<"npattern*coverage="<<totalCoverage<<"\n"; + // save sector index as TTree + return error; +} +#endif + +int FTKPatternBySectorWriter::AppendPatternsFromASCIIstream +(istream &in) { + // read patterns from a stream in ASCII format and append each pattern + // to the TTree, determined by the sector number + int error=0; + Long64_t nPattern; + int nLayerFile; + // first two words are + in>>nPattern>>nLayerFile; + if(in.fail()) { + Error("AppendPatternsFromASCIIstream") + <<"could not read number of patterns and layers\n"; + error++; + nPattern=0; + } else { + SetNLayers(nLayerFile); + Info("AppendPatternsFromASCIIstream") + <<"nPattern="<<nPattern<<" nLayer="<<GetNLayers()<<"\n"; + } + FTKPatternWithCoverage pattern(GetNLayers()); + int sector; + fDir.cd(); + int nread=0; + for(Long64_t i=0;i<nPattern;i++) { + Long64_t number; + Int_t data,coverage; + if(!(i%100000)) { + ShowProgress(TString::Format("%d00000",(int)(i/100000))); + } + in>>number; + for(int iLayer=0;iLayer<GetNLayers();iLayer++) { + in>>data; + pattern.GetHitPattern().SetHit(iLayer,data); + } + in>>sector>>coverage; + pattern.SetCoverage(coverage); + if(in.fail()) { + error++; + Warning("AppendPatternsFromASCIIstream") + <<"reading failed pos="<<i<<"\n"; + break; + } + error =AppendPattern(sector,pattern); + if(error) { + Error("AppendPatternsFromASCIIstream") + <<"writing to root file failed pos="<<i<<"\n"; + } + if(error) break; + nread++; + } + if(nread!=nPattern) { + Warning("AppendPatternsFromASCIIstream") + <<"only "<<nread<<"/"<<nPattern<<" patterns imported successfully\n"; + } + return error; +} + +int FTKPatternBySectorWriter::AppendMergedPatternsSector +(int sector,FTKPatternOneSectorOrdered const *orderedByCoverage) { + int error=0; + for(FTKPatternOneSectorOrdered::Ptr_t i=orderedByCoverage->Begin(); + i!=orderedByCoverage->End();i++) { + error=AppendPattern(sector,orderedByCoverage->GetCoverage(i), + orderedByCoverage->GetHitPattern(i)); + if(error) { + Error("AppendMergedPatternsSector")<<"error="<<error<<"\n"; + break; + } + } + return error; +} + +//================== class FTKPatternBySectorForestReader ==== -FTKPatternBySectorReader::FTKPatternBySectorReader(FTKRootFileChain &chain) - : FTKPatternBySectorBase("FTKPatternBySectorReader") { - int error(0); // not that this gets incremented but is currently not used! +FTKPatternBySectorForestReader::FTKPatternBySectorForestReader(FTKRootFileChain &chain,int *error) + : FTKPatternBySectorReader("FTKPatternBySectorForestReader") { fSourceName= TString::Format("FTKRootFileChain(length=%d",chain.GetLength()); for(int i=0;i<chain.GetLength();i++) { @@ -90,14 +493,17 @@ FTKPatternBySectorReader::FTKPatternBySectorReader(FTKRootFileChain &chain) for(int i=0;i<chain.GetLength();i++) { // open root file TString inputFileName=chain[i]; - TDirectory *in=OpenRootFileReadonly(inputFileName); + TDirectory *in= + FTKRootFile::Instance()->OpenRootFileReadonly(inputFileName); if(!in) { - Error("(constructor)") - <<"problem to open root file \""<<inputFileName<<"\", skipping\n"; - error++; + if(error) (*error)++; + else + Error("(constructor)") + <<"problem to open root file \""<<inputFileName<<"\", skipping\n"; } else { // read Sector information - FTKPatternBySectorReader *input=new FTKPatternBySectorReader(*in); + FTKPatternBySectorForestReader *input= + new FTKPatternBySectorForestReader(*in); SetNLayers(input->GetNLayers()); // store number of patterns per sector for(int sector=input->GetFirstSector();sector>=0; @@ -124,34 +530,13 @@ FTKPatternBySectorReader::FTKPatternBySectorReader(FTKRootFileChain &chain) delete in; } } + InitializeSectorByCoverageTable(); } -void FTKPatternBySectorReader::SelectSubregion(int nSub,int iSub) { - if(nSub>1) { - Info("SelectSubregion") - <<"splitting subregions nSub="<<nSub<<" iSub="<<iSub<<"\n"; - int nAll=fPatterns.size(); - for(PatternTreeBySector_t::iterator iSector=fPatterns.begin(); - iSector!=fPatterns.end();) { - int sector=(*iSector).first; - PatternTreeBySector_t::iterator i0=iSector; - iSector++; - if((sector%nSub)!=iSub) { - if((*i0).second) { - delete (*i0).second; - (*i0).second=0; - } - fPatterns.erase(i0); - } - } - Info("SelectSubregion")<<"selected "<<fPatterns.size()<<"/"<<nAll - <<" sectors\n"; - } -} - -FTKPatternBySectorReader::FTKPatternBySectorReader -(TDirectory &dir,const char *filename,std::set<int> const *sectorlist) - : FTKPatternBySectorBase(filename) { +FTKPatternBySectorForestReader::FTKPatternBySectorForestReader +(TDirectory &dir,set<int> const *sectorlist, + int *error) + : FTKPatternBySectorReader("FTKPatternBySectorForestReader") { // locate all TTree objects with pattern data // and store them in the STL map // fPatterns @@ -163,6 +548,8 @@ FTKPatternBySectorReader::FTKPatternBySectorReader while((o=next())) { TString name=((TObjString *)o)->GetString(); int sector=FTKPatternRootTreeReader::ExtractSectorNumber(name); + // for testing + //if(sector>100) continue; if((sector>=0)&&((!sectorlist)|| (sectorlist->find(sector)!=sectorlist->end()))) { FTKPatternRootTreeReader * &patternTree=fPatterns[sector]; @@ -181,14 +568,50 @@ FTKPatternBySectorReader::FTKPatternBySectorReader SetContentType(CONTENT_MERGED); } } + InitializeSectorByCoverageTable(); + if(GetContentType()==CONTENT_EMPTY) { + if(error) { + (*error)++; + } else { + Warning("FTKPatternBySectorForestReader") + <<"no pattern bank found in "<<dir.GetName()<<"\n"; + } + } } -FTKPatternBySectorReader::~FTKPatternBySectorReader() { +void FTKPatternBySectorForestReader::InitializeSectorByCoverageTable(void) { + if(GetContentType()!=CONTENT_NOTMERGED) { + int nSector=0; + int nPattern=0; + for(PatternTreeBySectorRO_t::const_iterator i=fPatterns.begin(); + i!=fPatterns.end();i++) { + int sector=(*i).first; + FTKPatternRootTreeReader *reader=(*i).second; + reader->ReadCoverageOnly(true); + reader->Rewind(); + while(reader->ReadNextPattern()) { + m_SectorByCoverage + [reader->GetPattern().GetCoverage()][sector]++; + nPattern++; + } + reader->Rewind(); + reader->ReadCoverageOnly(false); + nSector++; + } + Info("InitializeSectorByCoverageTable") + <<"nSector="<<nSector<<" nPattern="<<nPattern + <<" nCoverage="<<m_SectorByCoverage.size()<<"\n"; + } else { + Info("InitializeSectorByCoverageTable")<<"not ordered\n"; + } +} + +FTKPatternBySectorForestReader::~FTKPatternBySectorForestReader() { // clean up Long_t totalPatternsRead=0; int totalSectorsWithData=0; int totalSectors=0; - for(PatternTreeBySector_t::iterator i=fPatterns.begin(); + for(PatternTreeBySectorRO_t::iterator i=fPatterns.begin(); i!=fPatterns.end();i++) { Long_t nPatterns=(*i).second->GetNumberOfPatternsRead(); if(nPatterns>0) totalSectorsWithData++; @@ -205,194 +628,37 @@ FTKPatternBySectorReader::~FTKPatternBySectorReader() { } } -Long64_t FTKPatternBySectorReader::GetNPatterns(int sector) const { +Long64_t FTKPatternBySectorForestReader::GetNPatterns(int sector) const { // return number of patterns in one sector Long64_t n=0; - PatternTreeBySector_t::const_iterator i=fPatterns.find(sector); + PatternTreeBySectorRO_t::const_iterator i=fPatterns.find(sector); if(i != fPatterns.end()) { n=(*i).second->GetNPatterns(); } return n; } -Long64_t FTKPatternBySectorReader::GetNPatterns(void) const { - // return total number of patterns - Long64_t n=0; - for(PatternTreeBySector_t::const_iterator i=fPatterns.begin(); - i!=fPatterns.end();i++) { - n += (*i).second->GetNPatterns(); +int FTKPatternBySectorForestReader::GetFirstSector(void) const { + // return first (lowest) sector number + if(fPatterns.begin()!=fPatterns.end()) { + return (*fPatterns.begin()).first; + } else { + return -1; } - return n; } -void FTKPatternBySectorReader::GetNPatternsByCoverage -(std::map<int,int> &coverageMap) const { - if(GetContentType()==CONTENT_NOTMERGED) { - Fatal("GetNPatternsByCoverage")<<"patterns are not merged\n"; - } - //Info("GetNPatternsByCoverage")<<"number of sectors "<<fPatterns.size()<<"\n"; - for(PatternTreeBySector_t::const_iterator i=fPatterns.begin(); - i!=fPatterns.end();i++) { - FTKPatternRootTreeReader *reader=(*i).second; - reader->ReadCoverageOnly(true); - reader->Rewind(); - while(reader->ReadNextPattern()) { - coverageMap[reader->GetPattern().GetCoverage()]++; - } - reader->Rewind(); - reader->ReadCoverageOnly(false); - } -} - -void FTKPatternBySectorReader::GetNPatternsBySectorCoverage -(std::map<int,std::map<int,int> > §orCoverageMap) const { - if(GetContentType()==CONTENT_NOTMERGED) { - Fatal("GetNPatternsBySectorCoverage")<<"patterns are not merged\n"; - } - //Info("GetNPatternsBySectorCoverage")<<"number of sectors "<<fPatterns.size()<<"\n"; - for(PatternTreeBySector_t::const_iterator i=fPatterns.begin(); - i!=fPatterns.end();i++) { - FTKPatternRootTreeReader *reader=(*i).second; - map<int,int> &coverageMap(sectorCoverageMap[(*i).first]); - reader->ReadCoverageOnly(true); - reader->Rewind(); - while(reader->ReadNextPattern()) { - coverageMap[reader->GetPattern().GetCoverage()]++; - } - reader->ReadCoverageOnly(false); - } -} - -int FTKPatternBySectorReader::GetFirstSector(void) const { - // return first (lowest) sector number - if(fPatterns.begin()!=fPatterns.end()) { - return (*fPatterns.begin()).first; - } else { - return -1; - } -} - -int FTKPatternBySectorReader::GetNextSector(int sector) const { - // return next (higher) sector number - // if there is none, return -1 - PatternTreeBySector_t::const_iterator i=fPatterns.upper_bound(sector); - if(i!=fPatterns.end()) - return (*i).first; - return -1; -} - -FTKPatternOneSector *FTKPatternBySectorReader::MergePatterns -(int sector) const { - // load and merge patterns from one sector into memory - FTKPatternOneSector *r=new FTKPatternOneSector(sector); - PatternTreeBySector_t::const_iterator sectorPtr=fPatterns.find(sector); - if(sectorPtr!=fPatterns.end()) { - FTKPatternRootTreeReader *patternTree=(*sectorPtr).second; - patternTree->Rewind(); - while(patternTree->ReadNextPattern()) { - r->AddPattern(patternTree->GetPattern()); - } - } - return r; -} - - -Long64_t FTKPatternBySectorReader::GetPatternByCoverage(SectorSet_t& sectorByCoverage, int iSub, int nSub) { - //! Get patterns ordered by coverage - //! return total number of patterns - if(GetContentType()==CONTENT_NOTMERGED) { - Error("GetPatternByCoverage") - <<"data still have to be merged\n"; - } - Long64_t nPatTotal=0; - for(PatternTreeBySector_t::iterator i=fPatterns.begin(); - i!=fPatterns.end();i++) { - int sector=(*i).first; - // // skip all but the selected subregion - if((sector %nSub)!=iSub) continue; - FTKPatternRootTreeReader* const reader=(*i).second; - nPatTotal += reader->GetNPatterns(); - reader->Rewind(); - if(reader->ReadNextPattern()) - sectorByCoverage.insert(i); - } - Info("GetPatternByCoverage")<<"number of sectors: " - <<sectorByCoverage.size()<<" total pattern count: "<<nPatTotal<<" (iSub="<<iSub<<", nSub="<<nSub<<")\n"; - return nPatTotal; -} - - -int FTKPatternBySectorReader::WritePatternsToASCIIstream(std::ostream &out,int iSub,int nSub) { - // read all patterns from all sectors - // write data ordered by - // coverage - // then sector number - // then keep the order of the hits in the input TTree - int error=0; - if(GetContentType()==CONTENT_NOTMERGED) { - Error("WritePatternsToASCIIstream") - <<"data can not be saved as ASCII (still have to be merged)\n"; - } - - SectorSet_t sectorByCoverage; - Long64_t nPatTotal = GetPatternByCoverage(sectorByCoverage,iSub,nSub); - - // write header - out<<nPatTotal<<" "<<GetNLayers()<<"\n"; - - // write patterns - int number=0; - while(sectorByCoverage.begin()!=sectorByCoverage.end()) { - PatternTreeBySector_t::iterator oneSector=*sectorByCoverage.begin(); - sectorByCoverage.erase(sectorByCoverage.begin()); - FTKPatternRootTreeReader *tree=(*oneSector).second; - int sector=(*oneSector).first; - int coverage=tree->GetPattern().GetCoverage(); - //if(coverage<minCoverage) break; - bool isEmpty=false; - do { - if(!(number %1000)) { - ShowProgress(TString::Format("%d/%d/%d",number,sector,coverage)); - } - out<<number; - number++; - FTKHitPattern const &patternData=tree->GetPattern().GetHitPattern(); - for(unsigned iLayer=0;iLayer<patternData.GetSize();iLayer++) { - out<<' '<<patternData.GetHit(iLayer); - } - out<<' '<<sector<<' '<<coverage<<"\n"; - if(out.fail()) { - error++; - Error("WritePatternsToASCIIstream")<<"write operation failed\n"; - break; - } - if(!tree->ReadNextPattern()) { - isEmpty=true; - break; - } - } while(tree->GetPattern().GetCoverage()==coverage); - if(!isEmpty) { - sectorByCoverage.insert(oneSector); - } - } - Info("WritePatternsToASCIIstream") - <<" total patterns written: "<<number<<"\n"; - return error; -} - -//================== class FTKPatternBySectorBlockReader =================== - -FTKPatternBySectorBlockReader::FTKPatternBySectorBlockReader -(TDirectory &dir,std::set<int> const *sectorList) : - FTKPatternBySectorReader(dir,"FTKPatternBySectorBlockReader",sectorList) { -} - -void FTKPatternBySectorBlockReader::Rewind(void) { - if(GetContentType()==CONTENT_NOTMERGED) { - Fatal("Rewind")<<"patterns are not merged\n"; - } - for(PatternTreeBySector_t::iterator i=fPatterns.begin(); - i!=fPatterns.end();i++) { +int FTKPatternBySectorForestReader::GetNextSector(int sector) const { + // return next (higher) sector number + // if there is none, return -1 + PatternTreeBySectorRO_t::const_iterator i=fPatterns.upper_bound(sector); + if(i!=fPatterns.end()) + return (*i).first; + return -1; +} + +void FTKPatternBySectorForestReader::RewindSector(int sector) { + PatternTreeBySectorRO_t::iterator i=fPatterns.find(sector); + if(i!=fPatterns.end()) { FTKPatternRootTreeReader *reader=(*i).second; reader->Rewind(); reader->ReadNextPattern(); @@ -400,36 +666,38 @@ void FTKPatternBySectorBlockReader::Rewind(void) { } } -FTKPatternOneSector *FTKPatternBySectorBlockReader::Read +FTKPatternOneSector *FTKPatternBySectorForestReader::Read (int sector,int minCoverage) { - if(GetContentType()==CONTENT_NOTMERGED) { - Fatal("Rewind")<<"patterns are not merged\n"; + if((GetContentType()==CONTENT_NOTMERGED)&&(minCoverage>1)) { + Fatal("Read")<<"patterns are not merged\n"; } - FTKPatternOneSector *r=0; - PatternTreeBySector_t::iterator isector=fPatterns.find(sector); + FTKPatternOneSector *r=new FTKPatternOneSector(sector); + PatternTreeBySectorRO_t::iterator isector=fPatterns.find(sector); if(isector!=fPatterns.end()) { FTKPatternRootTreeReader *reader=(*isector).second; - do { + while(reader->ReadNextPattern()) { FTKPatternWithCoverage const &pattern=reader->GetPattern(); - if(pattern.GetCoverage()<minCoverage) break; - if(!r) { - r=new FTKPatternOneSector(sector); + if(pattern.GetCoverage()<minCoverage) { + Long64_t pos=reader->SeekBeg(-1); + reader->SeekBeg(pos-1); + break; } r->AddPattern(pattern); - } while( reader->ReadNextPattern()); + } reader->Suspend(); } return r; } -//================== class FTKPatternBySectorWriter ======================== +//================== class FTKPatternBySectorForestWriter ================== -FTKPatternBySectorWriter::FTKPatternBySectorWriter(TDirectory &dir) - : FTKPatternBySectorBase("FTKPatternBySectorWriter"),fDir(dir) { +FTKPatternBySectorForestWriter::FTKPatternBySectorForestWriter(TDirectory &dir) + : FTKPatternBySectorWriter("FTKPatternBySectorWriter",dir) { fSourceName=fDir.GetName(); + Info("FTKPatternBySectorForestWriter")<<"dir="<<fSourceName<<"\n"; } -FTKPatternBySectorWriter::~FTKPatternBySectorWriter() { +FTKPatternBySectorForestWriter::~FTKPatternBySectorForestWriter() { // clean up Info("~FTKPatternBySectorWriter")<<"saving trees\n"; int n=0; @@ -437,7 +705,7 @@ FTKPatternBySectorWriter::~FTKPatternBySectorWriter() { int totalSectorsWithData=0; int nPrint=0; TDirectory* lastDir=gDirectory; - for(PatternTreeBySector_t::iterator i=fPatterns.begin(); + for(PatternTreeBySectorRW_t::iterator i=fPatterns.begin(); i!=fPatterns.end();i++) { if(!(nPrint %100)) { ShowProgress(TString::Format("%d",(*i).first)); @@ -463,152 +731,7 @@ FTKPatternBySectorWriter::~FTKPatternBySectorWriter() { lastDir->cd(); // go back } -int FTKPatternBySectorWriter::AppendPatternsFromASCIIstream -(std::istream &in) { - // read patterns from a stream in ASCII format and append each pattern - // to the TTree, determined by the sector number - int error=0; - Long64_t nPattern; - int nLayerFile; - // first two words are - in>>nPattern>>nLayerFile; - if(in.fail()) { - Error("AppendPatternsFromASCIIstream") - <<"could not read number of patterns and layers\n"; - error++; - nPattern=0; - } else { - SetNLayers(nLayerFile); - Info("AppendPatternsFromASCIIstream") - <<"nPattern="<<nPattern<<" nLayer="<<GetNLayers()<<"\n"; - } - FTKPatternWithCoverage pattern(GetNLayers()); - int sector; - fDir.cd(); - int nread=0; - for(Long64_t i=0;i<nPattern;i++) { - Long64_t number; - Int_t data,coverage; - if(!(i%100000)) { - ShowProgress(TString::Format("%d00000",(int)(i/100000))); - } - in>>number; - for(int iLayer=0;iLayer<GetNLayers();iLayer++) { - in>>data; - pattern.GetHitPattern().SetHit(iLayer,data); - } - in>>sector>>coverage; - pattern.SetCoverage(coverage); - if(in.fail()) { - error++; - Warning("AppendPatternsFromASCIIstream") - <<"reading failed pos="<<i<<"\n"; - break; - } - error =AppendPattern(sector,pattern); - if(error) { - Error("AppendPatternsFromASCIIstream") - <<"writing to root file failed pos="<<i<<"\n"; - } - if(error) break; - nread++; - } - if(nread!=nPattern) { - Warning("AppendPatternsFromASCIIstream") - <<"only "<<nread<<"/"<<nPattern<<" patterns imported successfully\n"; - } - return error; -} - -int FTKPatternBySectorWriter::AppendMergedPatterns -(FTKPatternBySectorReader &source,int minCoverage) { - // loop over all sectors - // for each sector, import-and-merge patterns - int error=0; - SetNLayers(source.GetNLayers()); - // sector index, to be saved as extra TTree - map<int,map<int,Long64_t> > coverageSectorBeginEnd; - Info("AppendMergedPatterns")<<"loop over sectors\n"; - int nTotal=0; - int totalCoverage=0; - int printed=0; - double nextPrint=1.0; - for(int sector=source.GetFirstSector();sector>=0; - sector=source.GetNextSector(sector)) { - //ShowProgress(TString::Format("%d",sector)); - FTKPatternOneSector *mergedPatterns=source.MergePatterns(sector); - nTotal += mergedPatterns->GetNumberOfPatterns(); - FTKPatternOneSectorOrdered *orderedByCoverage= - mergedPatterns->OrderPatterns(FTKPatternOrderByCoverage(minCoverage)); - totalCoverage += orderedByCoverage->GetSummedCoverage(); - printed++; - if(printed>=nextPrint) { - ShowProgress(TString::Format("%d %d %d %d",sector, - mergedPatterns->GetNumberOfPatterns(), - nTotal,totalCoverage)); - nextPrint=printed*1.3; - } - error=AppendMergedPatterns(sector,orderedByCoverage); - if(error) {delete orderedByCoverage; delete mergedPatterns; break;} - int coverage=0; - int begin=0; - int nPat=0; - for(FTKPatternOneSectorOrdered::Ptr_t i=orderedByCoverage->Begin(); - i!=orderedByCoverage->End();i++) { - if(orderedByCoverage->GetCoverage(i)<minCoverage) break; - if(orderedByCoverage->GetCoverage(i)!=coverage) { - if(nPat!=begin) { - coverageSectorBeginEnd[coverage][sector]= - (((Long64_t)begin)<<32)|nPat; - } - begin=nPat; - coverage=orderedByCoverage->GetCoverage(i); - } - nPat++; - } - if(nPat!=begin) { - coverageSectorBeginEnd[coverage][sector]= - (((Long64_t)begin)<<32)|nPat; - - } - delete orderedByCoverage; - delete mergedPatterns; - } - Info(TString::Format("nsector=%d npattern=%d npattern*coverage=%d",printed, - nTotal,totalCoverage)); - // save sector index as TTree -#ifdef WRITE_INDEX_TREE - fDir.cd(); - TTree *tree=new TTree("sector_index","sector_index"); - int coverage,sector,begin,end; - tree->Branch("coverage",&coverage,"coverage/I"); - tree->Branch("sector",§or,"sector/I"); - tree->Branch("begin",&begin,"begin/I"); - tree->Branch("end",&end,"end/I"); - for(map<int,map<int,Long64_t> >::const_reverse_iterator - i=coverageSectorBeginEnd.rbegin();i!=coverageSectorBeginEnd.rend(); - i++) { - coverage=(*i).first; - int nPattern=0; - for(map<int,Long64_t>::const_iterator j=(*i).second.begin(); - j!=(*i).second.end();j++) { - sector=(*j).first; - begin=(*j).second>>32; - end=(*j).second; - nPattern += end-begin; - tree->Fill(); - } - Info("summary") - <<"coverage="<<coverage<<" nsector="<<(*i).second.size() - <<" npattern="<<nPattern<<"\n"; - } - tree->AutoSave("SaveSelf"); - delete tree; -#endif - return error; -} - -int FTKPatternBySectorWriter::AppendMergedPatterns +int FTKPatternBySectorForestWriter::AppendMergedPatternsSector (int sector,FTKPatternOneSectorOrdered const *orderedByCoverage) { // order patterns of one sector by coverage, then store them to disk int error=0; @@ -633,8 +756,8 @@ int FTKPatternBySectorWriter::AppendMergedPatterns return error; } -int FTKPatternBySectorWriter::AppendPattern -(int sector,FTKPatternWithCoverage const &pattern) { +int FTKPatternBySectorForestWriter::AppendPattern +(int sector,int coverage,FTKHitPattern const &pattern) { // add one pattern to the output int error=0; FTKPatternRootTree * &patternTree=fPatterns[sector]; @@ -646,7 +769,630 @@ int FTKPatternBySectorWriter::AppendPattern } else { SetContentType(CONTENT_NOTMERGED); } - error= !patternTree->WritePattern(pattern); + error= !patternTree->WritePattern(pattern,coverage); + return error; +} + +//================== class FTKPatternBySectorIndexed ================== + +/* int FTKPatternBySectorIndexed::GetBitsPerPlane(void) const { + return 64/GetNLayers(); + } */ + +char const *FTKPatternBySectorIndexed::s_cdataTreeName="CpatternData"; +char const *FTKPatternBySectorIndexed::s_indexTreeName="CpatternIndex"; + +char const *FTKPatternBySectorIndexed::s_cnpattBranchName="npatt"; +char const *FTKPatternBySectorIndexed::s_cnpattBranchID="npatt/i"; + +char const *FTKPatternBySectorIndexed::GetCPatternBranchName(uint32_t plane) { + while(s_cpatternBranch.size()<=plane) { + int n=s_cpatternBranch.size(); + s_cpatternBranch.push_back(Form("cpattern%d",n)); + } + return s_cpatternBranch[plane].c_str(); +} + +char const *FTKPatternBySectorIndexed::GetCPatternBranchID(uint32_t plane) { + while(s_cpatternID.size()<=plane) { + int n=s_cpatternID.size(); + s_cpatternID.push_back(Form("cpattern%d[npatt]/b",n)); + } + return s_cpatternID[plane].c_str(); +} + +char const *FTKPatternBySectorIndexed::s_sectorBranchName="sector"; +char const *FTKPatternBySectorIndexed::s_sectorBranchID="sector/i"; + +char const *FTKPatternBySectorIndexed::s_ncovBranchName="ncov"; +char const *FTKPatternBySectorIndexed::s_ncovBranchID="ncov/i"; + +char const *FTKPatternBySectorIndexed::s_npatternBranchName="npattern"; +char const *FTKPatternBySectorIndexed::s_npatternBranchID="npattern[ncov]/i"; + +char const *FTKPatternBySectorIndexed::s_coverageBranchName="coverage"; +char const *FTKPatternBySectorIndexed::s_coverageBranchID="coverage[ncov]/i"; + +vector<string> FTKPatternBySectorIndexed::s_cpatternBranch, + FTKPatternBySectorIndexed::s_cpatternID, + FTKPatternBySectorIndexed::s_ndecoderBranch, + FTKPatternBySectorIndexed::s_decoderDataBranch, + FTKPatternBySectorIndexed::s_ndecoderID, + FTKPatternBySectorIndexed::s_decoderDataID; + +char const *FTKPatternBySectorIndexed::GetNDecoderBranchName(uint32_t plane) { + while(s_ndecoderBranch.size()<=plane) { + int n=s_ndecoderBranch.size(); + s_ndecoderBranch.push_back(Form("ndecoder%d",n)); + } + return s_ndecoderBranch[plane].c_str(); +} + +char const *FTKPatternBySectorIndexed::GetNDecoderBranchID(uint32_t plane) { + while(s_ndecoderID.size()<=plane) { + int n=s_ndecoderID.size(); + s_ndecoderID.push_back(Form("ndecoder%d/i",n)); + } + return s_ndecoderID[plane].c_str(); +} + +char const *FTKPatternBySectorIndexed::GetDecoderDataBranchName(uint32_t plane) { + while(s_decoderDataBranch.size()<=plane) { + int n=s_decoderDataBranch.size(); + s_decoderDataBranch.push_back(Form("decoderData%d",n)); + } + return s_decoderDataBranch[plane].c_str(); +} + +char const *FTKPatternBySectorIndexed::GetDecoderDataBranchID(uint32_t plane) { + while(s_decoderDataID.size()<=plane) { + int n=s_decoderDataID.size(); + s_decoderDataID.push_back(Form("decoderData%d[ndecoder%d]/i",n,n)); + } + return s_decoderDataID[plane].c_str(); +} + +//================== class FTKPatternBySectorIndexedWriter ================== + +FTKPatternBySectorIndexedWriter::FTKPatternBySectorIndexedWriter +(TDirectory &dir) + : FTKPatternBySectorWriter("FTKPatternBySectorIndexedWriter",dir) { + fSourceName=fDir.GetName(); + Info("FTKPatternBySectorIndexedWriter")<<"dir="<<fSourceName<<"\n"; +m_cpatternDataTree=0; + m_cpatternIndexTree=0; + m_patternStorage=0; +} + +FTKPatternBySectorIndexedWriter::~FTKPatternBySectorIndexedWriter() { + Flush(); + fDir.cd(); + m_cpatternDataTree->Write(0,TObject::kOverwrite); + m_cpatternIndexTree->Write(0,TObject::kOverwrite); + delete m_cpatternDataTree; + delete m_cpatternIndexTree; +} + +int FTKPatternBySectorIndexedWriter::AppendPattern +(int sector,int coverage,FTKHitPattern const &pattern) { + // try two times + // at the end of each try, flush the table, so there is enough space + // for the second attempt + int error; + static int nPattern=0; + static int flush1=0; + static int flush2=0; + nPattern++; + for(int itry=0;itry<2;itry++) { + if(!m_patternStorage) + m_patternStorage=new CodedPatternStorage(GetNLayers()); + map<int,PatternData_t>::iterator iSector=m_patternData.find(sector); + if(iSector==m_patternData.end()) { + iSector=m_patternData.insert + (make_pair(sector,PatternData(GetNLayers(),m_patternStorage))) + .first; + } + PatternData_t §orData=(*iSector).second; + error=0; + uint32_t iCodedPattern=m_patternStorage->AllocatePattern(); + if(!iCodedPattern) { + // out of memory + error=2; + flush2++; + } + if(!error) { + // encode pattern for each layer and store as 64-bit word + for(int plane=GetNLayers()-1;plane>=0;plane--) { + // locate SSID in encoding table + int ssid=pattern.GetHit(plane); + map<int,int>::iterator iSSID= + sectorData.m_encoder[plane].second.find(ssid); + if(iSSID==sectorData.m_encoder[plane].second.end()) { + // not found, have to enlarge encoding table + // determine new coded ID + size_t id=sectorData.m_encoder[plane].second.size(); + if(id<sectorData.m_encoder[plane].first.size()) { + sectorData.m_encoder[plane].first[id]=ssid; + sectorData.m_encoder[plane].second[ssid]=id; + } else { + // encoding table full, have to flush + error=1; + flush1++; + break; + } + iSSID=sectorData.m_encoder[plane].second.find(ssid); + } + m_patternStorage->Store(iCodedPattern,plane,(*iSSID).second); + } + } + if(!error) { + CPatternWithCoverage_t::iterator iData= + sectorData.m_data.find(iCodedPattern); + if(iData==sectorData.m_data.end()) { + // insert new pattern + sectorData.m_data[iCodedPattern]+=coverage; + } else { + // increase coverage + (*iData).second+=coverage; + m_patternStorage->Delete(iCodedPattern); + } + } + if(error) { + // flush buffers + Flush(); + } else { + break; + } + } + /* static int loop=100; + if(loop) { + loop--; + if(!loop) exit(0); + } */ return error; } +void FTKPatternBySectorIndexedWriter::Flush(void) { + Info("Flush")<<"Number of sectors: "<<m_patternData.size()<<"\n"; + fDir.cd(); + if(!m_cpatternDataTree) { + + m_cpatternDataTree=new TTree(s_cdataTreeName,s_cdataTreeName); + m_cnpattBranch=m_cpatternDataTree->Branch + (s_cnpattBranchName,0,s_cnpattBranchID); + for(int plane=0;plane<GetNLayers();plane++) { + m_cpatternBranch.push_back(m_cpatternDataTree->Branch + (GetCPatternBranchName(plane),0,GetCPatternBranchID(plane))); + } + + m_cpatternIndexTree=new TTree(s_indexTreeName,s_indexTreeName); + m_cpatternIndexTree->Branch(s_sectorBranchName,&m_treeSector, + s_sectorBranchID); + m_cpatternIndexTree->Branch(s_ncovBranchName,&m_treeNCov, + s_ncovBranchID); + m_npatternBranch= + m_cpatternIndexTree->Branch(s_npatternBranchName,0, + s_npatternBranchID); + m_coverageBranch= + m_cpatternIndexTree->Branch(s_coverageBranchName,0, + s_coverageBranchID); + m_decoderDataBranch.resize(GetNLayers()); + m_decoderSize.resize(GetNLayers()); + for(int i=0;i<GetNLayers();i++) { + m_cpatternIndexTree->Branch + (GetNDecoderBranchName(i),&m_decoderSize[i], + GetNDecoderBranchID(i)); + m_decoderDataBranch[i]=m_cpatternIndexTree->Branch + (GetDecoderDataBranchName(i),0,GetDecoderDataBranchID(i)); + } + } + + // dump complete table from memory to disk + + for(map<int,PatternData_t>::iterator iSector= + m_patternData.begin();iSector!=m_patternData.end();iSector++) { + // dump table for one sector to disk + PatternData_t const & patternData=(*iSector).second; + m_treeSector=(*iSector).first; + + // order pattern data by coverage + map<uint32_t,list<int32_t> > ordered; + // ordered.first : coverage + // ordered.second : set(patterns) + for(CPatternWithCoverage_t::const_iterator i=patternData.m_data.begin(); + i!=patternData.m_data.end();i++) { + ordered[(*i).second].push_back((*i).first); + } + // + m_treeNCov=0; + int n=ordered.size(); + m_treeCoverage.resize(n); + m_treeNPattern.resize(n); + // reverse loop over coverage + for(map<uint32_t,list<int32_t> >::reverse_iterator iCov=ordered.rbegin(); + iCov!=ordered.rend();++iCov) { + // save data + m_treeCoverage[m_treeNCov]=(*iCov).first; + uint32_t &npatt=m_treeNPattern[m_treeNCov]; + npatt=(*iCov).second.size(); + m_cnpattBranch->SetAddress(&npatt); + vector<vector<uint8_t> > cpattern(GetNLayers()); + for(size_t plane=0;plane<cpattern.size();plane++) { + cpattern[plane].resize(npatt); + m_cpatternBranch[plane]->SetAddress(cpattern[plane].data()); + } + npatt=0; + for(list<int32_t>::const_iterator iPat=(*iCov).second.begin(); + iPat!=(*iCov).second.end();++iPat) { + for(size_t plane=0;plane<cpattern.size();plane++) { + cpattern[plane][npatt]=m_patternStorage->Get((*iPat),plane); + } + ++npatt; + } + m_cpatternDataTree->Fill(); + ++m_treeNCov; + } + // set branch addresses + m_npatternBranch->SetAddress(m_treeNPattern.data()); + m_coverageBranch->SetAddress(m_treeCoverage.data()); + // save encoder + for(int i=0;i<GetNLayers();i++) { + m_decoderDataBranch[i]->SetAddress + ((*iSector).second.m_encoder[i].first.data()); + m_decoderSize[i]=(*iSector).second.m_encoder[i].second.size(); + } + m_cpatternIndexTree->Fill(); + + } + // clear memory + m_patternData.clear(); + delete m_patternStorage; + m_patternStorage=0; +} + + + +//================== class FTKPatternBySectorIndexedReader ================== + +int FTKPatternBySectorIndexedReader::InitializeIndex +(TTree *indexTree,set<int> const *sectorList) { + int error=0; + // number of decoders + int nIndex=indexTree->GetEntries(); + Info("InitializeIndex")<<"Loading index, nIndex="<<nIndex<<"\n"; + // determine number of planes + uint32_t nPlane=0; + for(nPlane=0;indexTree->GetBranch(GetNDecoderBranchName(nPlane));++nPlane); + SetNLayers(nPlane); + m_patternData=new FTKPatternWithCoverage(nPlane); + + // one entry + typedef struct DecoderWithIndex { + // size of decoders per plane + uint32_t *m_decoderSize; + // number of Indexes + uint32_t m_ncov; + // sector number + uint32_t m_sectorID; + } DecoderWithIndex_t; + // full table, read from TTree + vector<DecoderWithIndex_t> decoderWithIndexTable(nIndex); + // set up pointers to decoder per plane + vector<uint32_t> decoderSizeData(nIndex*nPlane); + uint32_t *decoderSizePtr=decoderSizeData.data(); + for(int iIndex=0;iIndex<nIndex;++iIndex) { + decoderWithIndexTable[iIndex].m_decoderSize=decoderSizePtr; + decoderSizePtr += nPlane; + } + // read sector structure, decoder sizes and ncov + vector<uint32_t> ndecoder(nPlane); + uint32_t sector=666,ncov=666666; + for(uint32_t plane=0;plane<nPlane;++plane) { + indexTree->SetBranchAddress(GetNDecoderBranchName(plane), + &ndecoder[plane]); + } + indexTree->SetBranchAddress(s_sectorBranchName,§or); + indexTree->SetBranchAddress(s_ncovBranchName,&ncov); + + // first loop over the tree -> ncov,sector,ndecoder%d + for(int iIndex=0;iIndex<nIndex;iIndex++) { + indexTree->GetEntry(iIndex); + DecoderWithIndex_t &decoderWithIndex=decoderWithIndexTable[iIndex]; + decoderWithIndex.m_sectorID=sector; + decoderWithIndex.m_ncov=ncov; + for(uint32_t plane=0;plane<nPlane;++plane) { + decoderWithIndex.m_decoderSize[plane]=ndecoder[plane]; + } + } + uint32_t maxNCov=0; + uint32_t decoderSize=0; + for(int iIndex=0;iIndex<nIndex;iIndex++) { + DecoderWithIndex_t const &decoderWithIndex= + decoderWithIndexTable[iIndex]; + for(uint32_t plane=0;plane<nPlane;++plane) { + decoderSize += decoderWithIndex.m_decoderSize[plane]; + } + if(decoderWithIndex.m_ncov>maxNCov) maxNCov=decoderWithIndex.m_ncov; + } + Info("InitializeIndex") + <<"decoder size="<<decoderSize<<" maxNCov="<<maxNCov<<" \n"; + m_decoderData.resize(decoderSize); + m_decoderDataPtr.resize(nIndex*nPlane); + + vector<uint32_t> bufferCoverage(maxNCov); + vector<uint32_t> bufferNPattern(maxNCov); + // read index and decoder data and build up sector map + indexTree->SetBranchAddress(s_npatternBranchName,bufferNPattern.data()); + indexTree->SetBranchAddress(s_coverageBranchName,bufferCoverage.data()); + + uint32_t lastSector=0; + SetContentType(CONTENT_MERGED); + uint32_t cpatternEntry=0; + Info("InitializeIndex")<<"second loop\n"; + decoderSize=0; + for(int iIndex=0;iIndex<nIndex;iIndex++) { + DecoderWithIndex_t const &decoderWithIndex= + decoderWithIndexTable[iIndex]; + for(uint32_t plane=0;plane<nPlane;++plane) { + m_decoderDataPtr[iIndex*nPlane+plane]=&m_decoderData[decoderSize]; + indexTree->SetBranchAddress(GetDecoderDataBranchName(plane), + m_decoderDataPtr[iIndex*nPlane+plane]); + decoderSize+= decoderWithIndex.m_decoderSize[plane]; + } + indexTree->GetEntry(iIndex); + + uint32_t sector=decoderWithIndex.m_sectorID; + SectorInfo_t §orData=m_sectorTable[sector]; + if(!sectorData.m_coverageTable.size()) { + sectorData.m_nPattern=0; + } else if(lastSector!=sector) { + SetContentType(CONTENT_NOTMERGED); + } + // loop over coverage,npattern + for(uint32_t icov=0;icov<decoderWithIndex.m_ncov;++icov) { + DataBlock_t indexData; + indexData.m_decoderDataPtr=&m_decoderDataPtr[iIndex*nPlane]; + indexData.m_cpatternEntry=cpatternEntry; + indexData.m_nPattern=bufferNPattern[icov]; + ++cpatternEntry; + // skip sectors which are not in the sector list + if((!sectorList)||(sectorList->find(sector)!=sectorList->end())) { + sectorData.m_coverageTable.insert + (make_pair(bufferCoverage[icov],indexData)); + sectorData.m_coveragePtr=sectorData.m_coverageTable.rbegin(); + sectorData.m_nPattern+=indexData.m_nPattern; + } + } + lastSector=sector; + } + if(m_sectorTable.empty()) { + SetContentType(CONTENT_EMPTY); + error=1; + } + // if merged, set up coverage table + uint64_t nPattern=0; + for(SectorTable_t::const_iterator + iSector=m_sectorTable.begin(); + iSector!=m_sectorTable.end();iSector++) { + CoverageTable_t const &ctable=(*iSector).second.m_coverageTable; + for(CoverageTable_t::const_iterator iCov=ctable.begin(); + iCov!=ctable.end();iCov++) { + m_SectorByCoverage[(*iCov).first][(*iSector).first]+= + (*iCov).second.m_nPattern; + nPattern+=(*iCov).second.m_nPattern; + } + } + Info("InitializeIndex") + <<"nSector="<<m_sectorTable.size() + <<" nPattern="<<nPattern + <<" nCoverage="<<m_SectorByCoverage.size()<<"\n"; + return error; +} + +FTKPatternBySectorIndexedReader::FTKPatternBySectorIndexedReader +(TDirectory &dir,set<int> const *sectorList,int *error) + : FTKPatternBySectorReader("FTKPatternBySectorIndexedReader") { + TTree *indexTree; + dir.GetObject(s_indexTreeName,indexTree); + m_dataTree=0; + if(indexTree && !InitializeIndex(indexTree,sectorList)) { + dir.GetObject(s_cdataTreeName,m_dataTree); + } + if(indexTree) { + delete indexTree; + } + if(!m_dataTree) { + if(error) { + *error=1; + } else { + Error("FTKPatternBySectorIndexedReader")<<"not implemented\n"; + } + } +} + +FTKPatternBySectorIndexedReader::FTKPatternBySectorIndexedReader +(FTKRootFileChain &chain,int *error) + : FTKPatternBySectorReader("FTKPatternBySectorIndexedReader") { + TChain indexChain(s_indexTreeName); + int errorCount=0; + for(int i=0;i<chain.GetLength();i++) { + if(!indexChain.Add(chain[i].c_str(),0)) ++errorCount; + } + m_dataTree=0; + if((!errorCount) && !InitializeIndex(&indexChain,0)) { + TChain *dataChain=new TChain(s_cdataTreeName); + for(int i=0;i<chain.GetLength();i++) { + if(!dataChain->Add(chain[i].c_str(),0)) ++errorCount; + } + if(!errorCount) { + m_dataTree=dataChain; + } else { + delete dataChain; + } + } + if(!m_dataTree) { + if(error) { + *error=1; + } else { + Error("FTKPatternBySectorIndexedReader")<<"not implemented\n"; + } + } +} + +FTKPatternBySectorIndexedReader::~FTKPatternBySectorIndexedReader() { + if(m_dataTree) delete m_dataTree; +} + +Long64_t FTKPatternBySectorIndexedReader::GetNPatterns(int sector) const { + Long64_t n=0; + SectorTable_t::const_iterator i= + m_sectorTable.find(sector); + if(i!=m_sectorTable.end()) { + n=(*i).second.m_nPattern; + } + return n; +} + +int FTKPatternBySectorIndexedReader::GetFirstSector(void) const { + int sector=-1; + if(m_sectorTable.begin()!=m_sectorTable.end()) { + sector=(*m_sectorTable.begin()).first; + } + return sector; +} + +int FTKPatternBySectorIndexedReader::GetNextSector(int sector) const { + SectorTable_t::const_iterator i= + m_sectorTable.upper_bound(sector); + int r=-1; + if(i!=m_sectorTable.end()) { + r=(*i).first; + } + return r; +} + +void FTKPatternBySectorIndexedReader::RewindSector(int sector) { + SectorTable_t::iterator iSector= + m_sectorTable.find(sector); + if(iSector!=m_sectorTable.end()) { + SectorInfo_t &info=(*iSector).second; + info.m_coveragePtr=info.m_coverageTable.rbegin(); + } +} + +FTKPatternOneSector *FTKPatternBySectorIndexedReader::Read +(int sector,int minCoverage) { + FTKPatternOneSector *r=new FTKPatternOneSector(sector); + SectorTable_t::iterator iSector= + m_sectorTable.find(sector); + int nread=0; + if(iSector!=m_sectorTable.end()) { + SectorInfo_t &info=(*iSector).second; + uint32_t cnpatt,cnpattMax=0; + m_dataTree->SetBranchAddress(s_cnpattBranchName,&cnpatt); + for(CoverageTable_t::const_reverse_iterator icov=info.m_coveragePtr; + icov!=info.m_coverageTable.rend();++icov) { + cnpatt=(*icov).second.m_nPattern; + if(cnpattMax<cnpatt) cnpattMax=cnpatt; + } + vector<vector<uint8_t> > patternData(GetNLayers()); + for(size_t plane=0;plane<patternData.size();plane++) { + patternData[plane].resize(cnpattMax); + m_dataTree->SetBranchAddress(GetCPatternBranchName(plane), + patternData[plane].data()); + } + while(info.m_coveragePtr !=info.m_coverageTable.rend()) { + int coverage=(*info.m_coveragePtr).first; + if(coverage<minCoverage) break; + // read all patterns + DataBlock_t const &block=(*info.m_coveragePtr).second; + uint32_t entry=block.m_cpatternEntry; + m_dataTree->GetEntry(entry); + // loop over patterns + for(uint32_t ipattern=0;ipattern<cnpatt;++ipattern) { + // decode SSIDs + for(int plane=0;plane<GetNLayers();plane++) { + m_patternData->SetHit + (plane,block.m_decoderDataPtr[plane] + [patternData[plane][ipattern]]); + } + m_patternData->SetCoverage(coverage); + r->AddPattern(*m_patternData); + nread++; + } + info.m_coveragePtr++; + } + } + return r; +} + +int FTKPatternBySectorIndexedReader::ReadRaw +(int firstSector,std::list<FTKPatternOneSector *> §orList, + uint64_t maxPattern) { + // define a work-plan + uint64_t patternEstimate=0; + int readSector; + int nSector=0; + uint32_t cnpattMax=0; + map<uint32_t,pair<int,CoverageTable_t::const_reverse_iterator> > IOordered; + for(readSector=firstSector;readSector>=0; + readSector=GetNextSector(readSector)) { + SectorTable_t::iterator iSector=m_sectorTable.find(readSector); + if(iSector!=m_sectorTable.end()) { + SectorInfo_t &info=(*iSector).second; + while(info.m_coveragePtr !=info.m_coverageTable.rend()) { + DataBlock_t const &block=(*info.m_coveragePtr).second; + IOordered[block.m_cpatternEntry]= + make_pair(readSector,info.m_coveragePtr); + info.m_coveragePtr++; + patternEstimate += block.m_nPattern; + if(cnpattMax<block.m_nPattern) cnpattMax=block.m_nPattern; + } + nSector++; + if(patternEstimate>maxPattern) { + readSector=GetNextSector(readSector); + break; + } + } + } + Info("ReadRaw") + <<"reading "<<IOordered.size()<<" chunks "<<nSector<<" sectors\n"; + + uint32_t cnpatt; + m_dataTree->SetBranchAddress(s_cnpattBranchName,&cnpatt); + vector<vector<uint8_t> > patternData(GetNLayers()); + for(size_t plane=0;plane<patternData.size();plane++) { + patternData[plane].resize(cnpattMax); + m_dataTree->SetBranchAddress(GetCPatternBranchName(plane), + patternData[plane].data()); + } + map<int, FTKPatternOneSector *> sectorMap; + for(map<uint32_t,pair<int,CoverageTable_t::const_reverse_iterator> > + ::const_iterator iChunk=IOordered.begin();iChunk!=IOordered.end(); + iChunk++) { + int sector=(*iChunk).second.first; + FTKPatternOneSector *§orData=sectorMap[sector]; + if(!sectorData) sectorData=new FTKPatternOneSector(sector); + CoverageTable_t::const_reverse_iterator iCT=(*iChunk).second.second; + DataBlock_t const &block=(*iCT).second; + uint32_t entry=block.m_cpatternEntry; + m_dataTree->GetEntry(entry); + // loop over patterns + for(uint32_t ipattern=0;ipattern<cnpatt;++ipattern) { + // decode SSIDs + for(int plane=0;plane<GetNLayers();plane++) { + m_patternData->SetHit + (plane,block.m_decoderDataPtr[plane] + [patternData[plane][ipattern]]); + } + m_patternData->SetCoverage((*iCT).first); + sectorData->AddPattern(*m_patternData); + } + } + for(map<int,FTKPatternOneSector *>::const_iterator i=sectorMap.begin(); + i!=sectorMap.end();i++) { + sectorList.push_back((*i).second); + } + return readSector; +} + diff --git a/Trigger/TrigFTK/TrigFTKSim/src/FTK_AMBank.cxx b/Trigger/TrigFTK/TrigFTKSim/src/FTK_AMBank.cxx index 54b4eb982a19299d92a8578d126e9039a0a1ed81..0eb634799461df9bd6ba21858a33855d9b00757a 100644 --- a/Trigger/TrigFTK/TrigFTKSim/src/FTK_AMBank.cxx +++ b/Trigger/TrigFTK/TrigFTKSim/src/FTK_AMBank.cxx @@ -575,12 +575,12 @@ int FTK_AMBank::readROOTBankSectorOrdered(TFile* pattfile, int maxpatt) int iSub=0; int nSub=1; - FTKPatternBySectorReader* preader = new FTKPatternBySectorReader(*pattfile); + FTKPatternBySectorForestReader* preader = new FTKPatternBySectorForestReader(*pattfile); if(!preader) { cerr<<"Error. Cannot read root-file pattern from: "<<pattfile->GetName()<<endl; return -1; } - m_npatterns = preader->GetNPatterns(); + m_npatterns = preader->GetNPatternsTotal(); setNPlanes(preader->GetNLayers()); if (maxpatt<0 || m_npatterns<maxpatt) { @@ -601,54 +601,46 @@ int FTK_AMBank::readROOTBankSectorOrdered(TFile* pattfile, int maxpatt) readBankInit(); - FTKPatternBySectorReader::SectorSet_t sectorByCoverage; - preader->GetPatternByCoverage(sectorByCoverage,iSub,nSub); - // long npatt = preader->GetPatternByCoverage(sectorByCoverage,iSub,nSub); - // cout<< "Found # "<<npatt<<" patterns."<<endl; - // loop over sectors recursively, convert to 'FTKPattern' and write into root-tree int patternID = 0; cout << "Reading : [" << flush; - while(sectorByCoverage.begin()!=sectorByCoverage.end() && patternID<m_npatterns ) { - if (patternID%ipatt_step==0) cout << patternID/ipatt_step << flush; - - FTKPatternBySectorReader::PatternTreeBySector_t::iterator oneSector=*sectorByCoverage.begin(); - sectorByCoverage.erase(sectorByCoverage.begin()); // remove this sector temporarily - - FTKPatternRootTreeReader* tree=oneSector->second; - int sector=oneSector->first; - int coverage=tree->GetPattern().GetCoverage(); - bool isNotEmpty=false; - do { - // //showstats(m_patternID, m_npatterns); - // inputpattern->setPatternID(m_patternID); - // inputpattern->setSectorID(sector); - // inputpattern->setCoverage(coverage); - //cout<<patternID<<" "; //hier - const FTKHitPattern& patternData = tree->GetPattern().GetHitPattern(); - for( int iplane=0 ; iplane<m_nplanes ; iplane++ ) { - m_patterns[_SSPOS(patternID,iplane)] = patternData.GetHit(iplane); + preader->Rewind(); + for(int coverage=preader->GetMaxCoverage();coverage; + coverage=preader->GetNextCoverage(coverage)) { + if(patternID>=m_npatterns) break; + for(int sector=preader->GetCoverageFirstSector(coverage); + sector>=0;sector=preader->GetCoverageNextSector(coverage,sector)) { + if((nSub>1)&&(sector%nSub!=iSub)) continue; + if(patternID>=m_npatterns) break; + + if (patternID%ipatt_step==0) cout << patternID/ipatt_step << flush; + + FTKPatternOneSector *data=preader->Read(sector,coverage); + for(FTKPatternOneSector::Ptr_t i=data->Begin();i!=data->End();i++) { + if(patternID>=m_npatterns) break; + + // //showstats(m_patternID, m_npatterns); + // inputpattern->setPatternID(m_patternID); + // inputpattern->setSectorID(sector); + // inputpattern->setCoverage(coverage); + //cout<<patternID<<" "; //hier + const FTKHitPattern& patternData =data->GetHitPattern(i); + + for( int iplane=0 ; iplane<m_nplanes ; iplane++ ) { + m_patterns[_SSPOS(patternID,iplane)] = patternData.GetHit(iplane); //cout<<m_patterns[_SSPOS(patternID,iplane)]<<" "; //hier - } - m_patterns[_SSPOS(patternID,m_nplanes)] = sector; - m_patternCoverage[patternID] = coverage; - m_totalCoverage += coverage; - (*m_sectorCoverage)[m_patterns[_SSPOS(patternID,m_nplanes)]] += m_patternCoverage[patternID]; + } + m_patterns[_SSPOS(patternID,m_nplanes)] = sector; + m_patternCoverage[patternID] = coverage; + m_totalCoverage += coverage; + (*m_sectorCoverage)[m_patterns[_SSPOS(patternID,m_nplanes)]] += m_patternCoverage[patternID]; //cout<<sector<<" "<<coverage<<" "<<endl;; //hier - patternID++; // increment the global pattern ID - isNotEmpty = tree->ReadNextPattern(); - // if(!tree->ReadNextPattern()) { - // isEmpty=true; - // break; - // } - } - while(tree->GetPattern().GetCoverage()==coverage && isNotEmpty && patternID<m_npatterns); - // add sector again to (coveraged ordered) set - if(isNotEmpty) { - sectorByCoverage.insert(oneSector); - } - } // end while(sectors) + patternID++; // increment the global pattern ID + } // end for (patterns) + delete data; + } // end for(sector) + } // end for(coverage) cout << "]" << endl; cout<<patternID<<" patterns read successfully."<<endl; diff --git a/Trigger/TrigFTK/TrigFTKSim/src/FTK_CompressedAMBank.cxx b/Trigger/TrigFTK/TrigFTKSim/src/FTK_CompressedAMBank.cxx index 34e1566e6bad32bb7abb0a84d2717b5f9b2a844d..4229e068c82f64535461f1cf05bae32041ef94cb 100644 --- a/Trigger/TrigFTK/TrigFTKSim/src/FTK_CompressedAMBank.cxx +++ b/Trigger/TrigFTK/TrigFTKSim/src/FTK_CompressedAMBank.cxx @@ -2852,23 +2852,23 @@ int FTK_CompressedAMBank::readSectorOrderedBank(const char *name, int maxpatts, int nSub,int numSectorMax,int nlamb) { // create partitions based on the number of subregions std::list<Partition> partitions; - TDirectory *infile=FTKRootFile::Instance()->OpenRootFileReadonly(name); - if(infile) { + TDirectory *TSPfile=FTKRootFile::Instance()->OpenRootFileReadonly(name); + if(TSPfile) { std::set<int> subregion; - { - FTKPatternBySectorReader reader(*infile,"FTKPatternBySectorReader"); - reader.SelectSubregion(nSub,getSubID()); - for(int sector=reader.GetFirstSector();sector>=0; - sector=reader.GetNextSector(sector)) { - subregion.insert(sector); - } - int nSectorMax=-1; - if(numSectorMax>0) { - nSectorMax=nSub>0 ? numSectorMax/nSub : numSectorMax; - } - partitions.push_back(Partition(maxpatts,nSectorMax,subregion)); + FTKPatternBySectorReader *TSPreader= + FTKPatternBySectorReader::Factory(*TSPfile); + for(int sector=TSPreader->GetFirstSector();sector>=0; + sector=TSPreader->GetNextSector(sector)) { + if((nSub>1)&&(sector%nSub!=getSubID())) continue; + subregion.insert(sector); } - delete infile; + delete TSPreader; + int nSectorMax=-1; + if(numSectorMax>0) { + nSectorMax=nSub>0 ? numSectorMax/nSub : numSectorMax; + } + partitions.push_back(Partition(maxpatts,nSectorMax,subregion)); + delete TSPfile; } else { Error("readSectorOrderedBank")<<"failed to open \""<<name<<"\"\n"; } @@ -3031,9 +3031,9 @@ int FTK_CompressedAMBank::readPartitionedSectorOrderedBank int nlamb) { // open file int error=0; - TDirectory *infile=FTKRootFile::Instance()->OpenRootFileReadonly(name); + TDirectory *TSPfile=FTKRootFile::Instance()->OpenRootFileReadonly(name); // read sector-ordered root file (TSP patterns) - if(infile) { + if(TSPfile) { int partitionCount=0; int sectorCount=0; int nLayer=0; @@ -3062,13 +3062,18 @@ int FTK_CompressedAMBank::readPartitionedSectorOrderedBank // // if the number of subregions is specified, it will ignore all // sectors which do not belong to the active subregion - FTKPatternBySectorBlockReader reader(*infile,&(*iPartition).fSectorSet); + FTKPatternBySectorReader *TSPreader=FTKPatternBySectorReader::Factory + (*TSPfile,&(*iPartition).fSectorSet); + if(TSPreader->GetNLayers()<=0) { + Fatal("readSectorOrderedBank") + <<"number of layers (reader) "<<TSPreader->GetNLayers()<<"\n"; + } if(!nLayer) { - setNPlanes(reader.GetNLayers()); - nLayer=reader.GetNLayers(); + setNPlanes(TSPreader->GetNLayers()); + nLayer=TSPreader->GetNLayers(); } Info("readSectorOrderedBank") - <<"reading from "<<name<<" nLayer="<<reader.GetNLayers() + <<"reading from "<<name<<" nLayer="<<TSPreader->GetNLayers() <<" partition "<<partitionCount<<" maxpattern=" <<(*iPartition).fNumPatternMax<<" maxSector="<<(*iPartition).fNumSectorMax <<"\n"; @@ -3077,10 +3082,10 @@ int FTK_CompressedAMBank::readPartitionedSectorOrderedBank // can be read printVmemUsage("before GetNPatternsByCoverage"); std::map<int,int> coverageMap; - reader.GetNPatternsByCoverage(coverageMap); + TSPreader->GetNPatternsByCoverage(coverageMap); printVmemUsage("after GetNPatternsByCoverage"); std::map<int,int>::const_reverse_iterator i=coverageMap.rbegin(); - reader.Rewind(); + TSPreader->Rewind(); printVmemUsage("after Rewind"); // print bank statistics // determine total number of patterns , patterns*coverage, <coverage> @@ -3117,8 +3122,8 @@ int FTK_CompressedAMBank::readPartitionedSectorOrderedBank // read all patterns in the given coverage range // (not reading the smallest coverage) // these patterns for sure will fit into the DC bank - //for(int sector=reader.GetFirstSector();sector>=0; - // sector=reader.GetNextSector(sector)) { + //for(int sector=TSPreader->GetFirstSector();sector>=0; + // sector=TSPreader->GetNextSector(sector)) { for(std::set<int>::const_iterator sectorPtr= (*iPartition).fSectorSet.begin(); sectorPtr!=(*iPartition).fSectorSet.end();sectorPtr++) { @@ -3127,7 +3132,7 @@ int FTK_CompressedAMBank::readPartitionedSectorOrderedBank if((int)loadedSectorList.size()>=maxNumSector) continue; loadedSectorList.insert(sector); } - FTKPatternOneSector *patterns=reader.Read(sector,covEnd+1); + FTKPatternOneSector *patterns=TSPreader->Read(sector,covEnd+1); if(patterns) { insertPatterns(sector,patterns,maxpatts,dcPatterns,nDC,nTSP); delete patterns; @@ -3143,7 +3148,7 @@ int FTK_CompressedAMBank::readPartitionedSectorOrderedBank if((int)loadedSectorList.size()>=maxNumSector) continue; loadedSectorList.insert(sector); } - FTKPatternOneSector *patterns=reader.Read(sector,covEnd); + FTKPatternOneSector *patterns=TSPreader->Read(sector,covEnd); if(patterns) { insertPatterns(sector,patterns,maxpatts,dcPatterns,nDC,nTSP); delete patterns; @@ -3162,6 +3167,7 @@ int FTK_CompressedAMBank::readPartitionedSectorOrderedBank <<"\n"; partitionCount++; sectorCount += loadedSectorList.size(); + delete TSPreader; } Info("readSectorOrderedBank") <<"total number of patterns read "<<nTSP<<" stored "<<nDC @@ -3211,13 +3217,11 @@ int FTK_CompressedAMBank::readPartitionedSectorOrderedBank delete [] ssidData; readBankPostprocessing("readSectorOrderedBank"); + delete TSPfile; } else { Error("readSectorOrderedBank")<<"failed to open \""<<name<<"\"\n"; error++; } - if(infile) { - delete infile; - } return error; } diff --git a/Trigger/TrigFTK/TrigFTKSim/src/tsp/TSPROOTBankGenerator.cxx b/Trigger/TrigFTK/TrigFTKSim/src/tsp/TSPROOTBankGenerator.cxx index a9c17d6ea8eb9c2a1f82213b813dbedc060e07db..59d47ca9a3bf51e1cb35711104e27b94219d26b7 100644 --- a/Trigger/TrigFTK/TrigFTKSim/src/tsp/TSPROOTBankGenerator.cxx +++ b/Trigger/TrigFTK/TrigFTKSim/src/tsp/TSPROOTBankGenerator.cxx @@ -63,12 +63,12 @@ TSPROOTBankGenerator::TSPROOTBankGenerator(FTKSetup* setup, const std::vector<FT FTKSetup::PrintMessageFmt(ftk::sevr, "Cannot read root-file from: %s\n", inputBank.c_str()); throw; } - m_preader = new FTKPatternBySectorReader(*file); + m_preader = new FTKPatternBySectorForestReader(*file); if(!m_preader) { FTKSetup::PrintMessageFmt(ftk::sevr, "Cannot read root-file pattern from: %s\n", inputBank.c_str()); throw; } - m_npatterns = m_preader->GetNPatterns(); + m_npatterns = m_preader->GetNPatternsTotal(); m_nplanes = m_preader->GetNLayers(); } else { @@ -114,30 +114,30 @@ void TSPROOTBankGenerator::generate() throw (TSPPatternReadException){ // the first step is just the conversion of the original (ASCII or root) bank a TTree if ( m_preader ) { - FTKPatternBySectorReader::SectorSet_t sectorByCoverage; - long npatt = m_preader->GetPatternByCoverage(sectorByCoverage,m_iSub,m_nSub); + long npatt = m_preader->GetNPatternsTotal(); // update number of patterns (for the current subregion) m_npatterns = m_maxPatterns >= 0 ? m_maxPatterns : npatt; // loop over sectors recursively, convert to 'FTKPattern' and write into root-tree - while(sectorByCoverage.begin()!=sectorByCoverage.end() && m_patternID<m_npatterns ) { - FTKPatternBySectorReader::PatternTreeBySector_t::iterator oneSector=*sectorByCoverage.begin(); - sectorByCoverage.erase(sectorByCoverage.begin()); // remove this sector temporarily - - FTKPatternRootTreeReader* tree=oneSector->second; - int sector=oneSector->first; - - int coverage=tree->GetPattern().GetCoverage(); - if ( coverage < m_mincoverage ) break; - - bool isNotEmpty=false; - do { + // coverage loop + for(int coverage=m_preader->GetMaxCoverage();coverage; + coverage=m_preader->GetNextCoverage(coverage)) { + if(m_patternID>=m_npatterns) break; + // sector loop + for(int sector=m_preader->GetCoverageFirstSector(coverage);sector>=0; + sector=m_preader->GetCoverageNextSector(coverage,sector)) { + if(m_patternID>=m_npatterns) break; + if ( coverage < m_mincoverage ) break; + FTKPatternOneSector *data=m_preader->Read(sector,coverage); + // pattern loop + for(FTKPatternOneSector::Ptr_t i=data->Begin();i!=data->End();i++) { + if(m_patternID>=m_npatterns) break; showstats(m_patternID, m_npatterns); inputpattern->setPatternID(m_patternID); inputpattern->setSectorID(sector); inputpattern->setCoverage(coverage); - const FTKHitPattern& patternData=tree->GetPattern().GetHitPattern(); + const FTKHitPattern& patternData=data->GetHitPattern(i); for(unsigned int iLayer=0;iLayer<patternData.GetSize();iLayer++) { inputpattern->setSSID(iLayer,patternData.GetHit(iLayer)); } @@ -146,22 +146,18 @@ void TSPROOTBankGenerator::generate() throw (TSPPatternReadException){ bank0->Fill(); m_patternID++; // increment the global pattern ID - isNotEmpty = tree->ReadNextPattern(); // if(!tree->ReadNextPattern()) { // isEmpty=true; // break; // } - } - while(isNotEmpty && tree->GetPattern().GetCoverage()==coverage && m_patternID<m_npatterns); - // add sector again to (coveraged ordered) set - if(isNotEmpty) { - sectorByCoverage.insert(oneSector); - } - } // end while(sectors) - - delete m_preader; - m_preader = NULL; - m_npatterns = m_patternID; + } // end loop over patterns + delete data; + } // end loop over sectors + } // end loop over coverage + + delete m_preader; + m_preader = NULL; + m_npatterns = m_patternID; } else { for (int ipatt = 0; ipatt < m_npatterns; ++ipatt) {