WIP: Algorithms for VXDTestbeam

Reconstruction/CMakeLists.txt

@@ -4,3 +4,4 @@ add_subdirectory(PFA)
 add_subdirectory(SiliconTracking)
 add_subdirectory(Tracking)
 add_subdirectory(RecGenfitAlg)
+add_subdirectory(VXDTestbeam)

Reconstruction/VXDTestbeam/CMakeLists.txt

@@ -0,0 +1,28 @@
+# Headers and Libraries
+
+# Modules
+gaudi_add_module(VXDTestbeam
+                 SOURCES src/TBDataReader/TBDataReader.cpp
+                 SOURCES src/TBDataReader/TBDataClusterReader.cpp
+                 LINK #DetInterface
+                      k4FWCore::k4FWCore
+                      Gaudi::GaudiAlgLib Gaudi::GaudiKernel
+                      ${LCIO_LIBRARIES} 
+                      ${DD4hep_COMPONENT_LIBRARIES}
+                      EDM4HEP::edm4hep EDM4HEP::edm4hepDict
+                      ${podio_LIBRARIES} podio::podioRootIO
+)
+
+if (GenFit_FOUND)
+  target_link_libraries(VXDTestbeam PUBLIC GenFit::genfit2)
+endif()
+
+target_include_directories(VXDTestbeam PUBLIC
+  $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}>/include
+  $<INSTALL_INTERFACE:${CMAKE_INSTALL_INCLUDEDIR}>)
+
+install(TARGETS VXDTestbeam
+  EXPORT CEPCSWTargets
+  RUNTIME DESTINATION "${CMAKE_INSTALL_BINDIR}" COMPONENT bin
+  LIBRARY DESTINATION "${CMAKE_INSTALL_LIBDIR}" COMPONENT shlib
+  COMPONENT dev)

Reconstruction/VXDTestbeam/options/TBDataClusterReader.py

@@ -0,0 +1,24 @@
+#!/usr/bin/env python
+
+from Gaudi.Configuration import *
+
+from Configurables import k4DataSvc
+dsvc = k4DataSvc("EventDataSvc")
+
+from Configurables import TBDataClusterReader
+alg = TBDataClusterReader("TBDataClusterReader")
+alg.TrackerHitOut.Path = "TrackerHit"
+
+from Configurables import PodioOutput
+out = PodioOutput("out")
+out.filename = "Testbeamcluster.root"
+out.outputCommands = ["keep *"]
+
+# ApplicationMgr
+from Configurables import ApplicationMgr
+ApplicationMgr( TopAlg = [alg, out],
+                EvtSel = 'NONE',
+                EvtMax = 1883792,
+                ExtSvc=[dsvc],
+                # OutputLevel=DEBUG
+)

Reconstruction/VXDTestbeam/options/TBDataReader.py

@@ -0,0 +1,24 @@
+#!/usr/bin/env python
+
+from Gaudi.Configuration import *
+
+from Configurables import k4DataSvc
+dsvc = k4DataSvc("EventDataSvc")
+
+from Configurables import TBDataReader
+alg = TBDataReader("TBDataReader")
+alg.RawTimeSeriesOut.Path = "RawTimeSeries"
+
+from Configurables import PodioOutput
+out = PodioOutput("out")
+out.filename = "Testbeamdata.root"
+out.outputCommands = ["keep *"]
+
+# ApplicationMgr
+from Configurables import ApplicationMgr
+ApplicationMgr( TopAlg = [alg, out],
+                EvtSel = 'NONE',
+                EvtMax = 8285,
+                ExtSvc=[dsvc],
+                # OutputLevel=DEBUG
+)

Reconstruction/VXDTestbeam/src/TBDataAnalysis/pixel.cpp

@@ -0,0 +1,71 @@
+// read csv
+#include <iostream>
+#include <fstream>
+#include <cstdlib>
+#include <sstream>
+
+// dependence
+#include "pixel.h"
+#include "edm4hep/RawTimeSeriesCollection.h"
+
+// root
+#include <TFile.h> 
+#include <TTree.h>
+
+DECLARE_COMPONENT(PixelWrite)
+
+PixelWrite::PixelWrite(const std::string& name, ISvcLocator* svcLoc)
+    : GaudiAlgorithm(name, svcLoc)
+{
+    declareProperty("RawTimeSeriesOut", m_hitCol);
+}
+
+StatusCode PixelWrite::initialize()
+{
+    info() << "begin initialize PixelWrite" << endmsg;
+
+    // get file path
+    tree_file = new TFile("/afs/ihep.ac.cn/users/l/lishuqi/public/TB/AllSameTimeHitsLoose2_9.root", "READ");
+
+    if (!tree_file->IsOpen()) {
+        info() << "failed open file" << endmsg;
+    }
+    tree = (TTree*)tree_file->Get("TimeInfo");
+
+    tree->SetBranchAddress("timeChip", &timeChip);
+    tree->SetBranchAddress("planeID", &planeID);
+    tree->SetBranchAddress("row", &row);
+    tree->SetBranchAddress("col", &col);
+
+    info() << "Total Number of Events = " << tree->GetEntries() << endmsg;
+    return GaudiAlgorithm::initialize();
+}
+
+StatusCode PixelWrite::execute()
+{
+    info() << "Executing Event " << eventid << endmsg;
+
+    tree->GetEntry(eventid);
+    auto hits = m_hitCol.createAndPut();
+    hits->setID(eventid);
+    int event_size = timeChip->size();
+    for (int i=0; i < event_size; i++){
+        auto hit = hits->create();
+        hit.setTime(timeChip->at(i));
+        // planeID:8,row:12,col:12
+        dd4hep::rec::CellID ncellid;
+        m_decoder.set(ncellid, "planeID", planeID->at(i));
+        m_decoder.set(ncellid, "row", row->at(i));
+        m_decoder.set(ncellid, "col", col->at(i));
+        hit.setCellID(ncellid);
+    }
+
+    eventid += 1;
+    return StatusCode::SUCCESS;
+}
+
+StatusCode PixelWrite::finalize()
+{
+    info() << "finalize PixelWrite" << endmsg;
+    return GaudiAlgorithm::finalize();
+}

Reconstruction/VXDTestbeam/src/TBDataAnalysis/pixel.h

@@ -0,0 +1,49 @@
+#ifndef PIXEL_H
+#define PIXEL_H
+
+#include "k4FWCore/DataHandle.h"
+#include "GaudiAlg/GaudiAlgorithm.h"
+#include "DD4hep/Detector.h"
+#include <DDRec/DetectorData.h>
+#include <DDRec/CellIDPositionConverter.h> 
+
+namespace edm4hep {
+    class RawTimeSeriesCollection;
+}
+
+namespace dd4hep {
+    namespace rec {
+        class CellIDPositionConverter;
+    }
+    namespace DDSegmentation {
+        class BitFieldCoder;
+    }
+}
+
+class PixelWrite : public GaudiAlgorithm
+{
+
+    public :
+
+        PixelWrite(const std::string& name, ISvcLocator* svcLoc);
+
+        virtual StatusCode initialize();
+        virtual StatusCode execute();
+        virtual StatusCode finalize();
+
+    private :
+
+        DataHandle<edm4hep::RawTimeSeriesCollection> m_hitCol{"RawTimeSeries", Gaudi::DataHandle::Writer, this};
+        dd4hep::DDSegmentation::BitFieldCoder m_decoder{"planeID:8,row:12,col:12"};
+        // dd4hep::DDSegmentation::BitFieldCoder* m_decoder;
+
+        int eventid = 0;
+        TFile *tree_file = NULL;
+        TTree *tree = NULL;
+        std::vector<int>* timeChip = NULL;
+        std::vector<int>* planeID = NULL;
+        std::vector<int>* row = NULL;
+        std::vector<int>* col = NULL;
+};
+
+#endif  // PIXEL_H

Reconstruction/VXDTestbeam/src/TBDataReader/TBDataClusterReader.cpp

@@ -0,0 +1,76 @@
+// read csv
+#include <iostream>
+#include <fstream>
+#include <cstdlib>
+#include <sstream>
+
+// dependence
+#include "TBDataClusterReader.h"
+#include "edm4hep/TrackerHitCollection.h"
+
+// root
+#include <TFile.h> 
+#include <TTree.h>
+
+DECLARE_COMPONENT(TBDataClusterReader)
+
+TBDataClusterReader::TBDataClusterReader(const std::string& name, ISvcLocator* svcLoc)
+    : GaudiAlgorithm(name, svcLoc)
+{
+    declareProperty("TrackerHitOut", m_hitCol);
+}
+
+StatusCode TBDataClusterReader::initialize()
+{
+    info() << "begin initialize TBDataClusterReader" << endmsg;
+
+    // get file path
+    tree_file = new TFile("/afs/ihep.ac.cn/users/l/lishuqi/scratchfs_2Weeks/pub/Clusterxy_all2_9.root", "READ");
+
+    if (!tree_file->IsOpen()) {
+        info() << "failed open file" << endmsg;
+    }
+    tree = (TTree*)tree_file->Get("Event");
+
+    tree->SetBranchAddress("planeID", &planeID);
+    tree->SetBranchAddress("col", &col);
+    tree->SetBranchAddress("row", &row);
+    tree->SetBranchAddress("hitsNum", &hitsNum);
+    tree->SetBranchAddress("hitsNumX", &hitsNumX);
+    tree->SetBranchAddress("hitsNumY", &hitsNumY);
+
+    info() << "Total Number of Events = " << tree->GetEntries() << endmsg;
+    return GaudiAlgorithm::initialize();
+}
+
+StatusCode TBDataClusterReader::execute()
+{
+    info() << "Executing Event " << eventid << endmsg;
+
+    tree->GetEntry(eventid);
+    auto hits = m_hitCol.createAndPut();
+    hits->setID(eventid);
+    int event_size = planeID->size();
+    for (int i=0; i < event_size; i++){
+        auto hit = hits->create();
+        // use float type storage
+        hit.setTime(col->at(i));
+        hit.setEDep(row->at(i));
+        // planeID:8,hitsNum:8,hitsNumX:8,hitsNumY:8
+        dd4hep::rec::CellID ncellid;
+        m_decoder.set(ncellid, "planeID", planeID->at(i));
+        m_decoder.set(ncellid, "hitsNum", hitsNum->at(i));
+        m_decoder.set(ncellid, "hitsNumX", hitsNumX->at(i));
+        m_decoder.set(ncellid, "hitsNumY", hitsNumY->at(i));
+        hit.setCellID(ncellid);
+    }
+
+    eventid += 1;
+    return StatusCode::SUCCESS;
+}
+
+StatusCode TBDataClusterReader::finalize()
+{
+    info() << "finalize TBDataClusterReader" << endmsg;
+    return GaudiAlgorithm::finalize();
+}

Reconstruction/VXDTestbeam/src/TBDataReader/TBDataClusterReader.h

@@ -0,0 +1,51 @@
+#ifndef TBDATACLUSTERREADER_H
+#define TBDATACLUSTERREADER_H
+
+#include "k4FWCore/DataHandle.h"
+#include "GaudiAlg/GaudiAlgorithm.h"
+#include "DD4hep/Detector.h"
+#include <DDRec/DetectorData.h>
+#include <DDRec/CellIDPositionConverter.h> 
+
+namespace edm4hep {
+    class TrackerHitCollection;
+}
+
+namespace dd4hep {
+    namespace rec {
+        class CellIDPositionConverter;
+    }
+    namespace DDSegmentation {
+        class BitFieldCoder;
+    }
+}
+
+class TBDataClusterReader : public GaudiAlgorithm
+{
+
+    public :
+
+        TBDataClusterReader(const std::string& name, ISvcLocator* svcLoc);
+
+        virtual StatusCode initialize();
+        virtual StatusCode execute();
+        virtual StatusCode finalize();
+
+    private :
+
+        DataHandle<edm4hep::TrackerHitCollection> m_hitCol{"TrackerHit", Gaudi::DataHandle::Writer, this};
+        dd4hep::DDSegmentation::BitFieldCoder m_decoder{"planeID:8,hitsNum:8,hitsNumX:8,hitsNumY:8"};
+        // dd4hep::DDSegmentation::BitFieldCoder* m_decoder;
+
+        int eventid = 0;
+        TFile *tree_file = NULL;
+        TTree *tree = NULL;
+        std::vector<int>* planeID = NULL;
+        std::vector<int>* hitsNum = NULL;
+        std::vector<int>* hitsNumX = NULL;
+        std::vector<int>* hitsNumY = NULL;
+        std::vector<float>* row = NULL;
+        std::vector<float>* col = NULL;
+};
+
+#endif  // TBDATACLUSTERREADER_H