CoMPASS_BinReader/Bin2Root2.cpp

/********************************************
 * 
 * This Bin2Root2.cpp is specially designed for CAEN CoMPASS output
 * 
 * that all boards and channels are combined together.
 * 
 * The data file name has format DataR_run_XXX_YYY.BIN
 * where XXX is run number, YYY from 1 when data size > 1 GB
 * The first file of the run is DataR_run_XXX.BIN
 * 
 * In this case, only the first .BIN file has header. 
 * 
 * Also, the timestamp is not fully sorted. 
 * 
 * Timeoffset is hardcoded.
 * 
 * Created by Ryan Tang, goluckyryan@gmail.com
 *                   last edit @ 2024-09-25
 * 
 *********************************************/

#include "BinReader.h"
#include <sys/time.h> /** struct timeval, select() */

inline unsigned int getTime_us(){
  unsigned int time_us;
  struct timeval t1;
  struct timezone tz;
  gettimeofday(&t1, &tz);
  time_us = (t1.tv_sec) * 1000 * 1000 + t1.tv_usec;
  return time_us;
}

#include "TFile.h"
#include "TTree.h"
#include "TMacro.h"

#define MAX_MULTI  100
#define debug 0 // for > 1 number of hit to debug

#include <sstream>
#include <string>
#include <vector>
#include <utility>
#include <filesystem>
#include <regex>

std::vector<std::pair<int, std::string>> findDataFiles(const std::string& patternPrefix) {
  std::vector<std::pair<int, std::string>> results;
  std::regex pattern1(patternPrefix + "_" + R"_((\d{1,3})\.BIN)_");
  std::regex pattern2(patternPrefix + R"(\.BIN)");

  for (const auto& entry : std::filesystem::directory_iterator(".")) {
    if (entry.is_regular_file()) {
      std::string filename = entry.path().filename().string();

      printf("%s\n", filename.c_str());
      std::smatch match;

      // Check for pattern DataR_run_34_YYY.BIN
      if (std::regex_match(filename, match, pattern1)) {
        int yyy = std::stoi(match[1].str());
        results.emplace_back(yyy, filename);
      } 
      // Check for pattern DataR_run_34.BIN
      else if (std::regex_match(filename, pattern2)) {
        results.emplace_back(0, filename);
      }
    }
  }

  std::sort(results.begin(), results.end(), [](const std::pair<int, std::string>& a, const std::pair<int, std::string>& b) {
    return a.first < b.first;
  });

  return results;
}

//^#############################################################
//^#############################################################
int main(int argc, char **argv) {
  
  printf("=========================================\n");
  printf("===      CoMPASS Binary to root       ===\n");
  printf("=========================================\n");  
  if (argc < 3)    {
    printf("Incorrect number of arguments:\n");
    printf("%s [timeWindow]  [RunNum] \n", argv[0]);
    printf("       timeWindow : in ns, -1 = no event building \n");   
    printf(" time Offset File : 0 for nothing \n");   
    printf("\n\n");

    return 1;
  }

  unsigned int runStartTime = getTime_us();

  ///============= read input
  long timeWindow = atoi(argv[1]);
  short runNum = atoi(argv[2]);

  ///============ time offset
  std::vector<int64_t> timeOffsetList =  {0, 0, 0, 0};

  ///=========== Form outFileName;
  TString outFileName = Form("run%03d_%s.root", runNum, (timeWindow < 0 ? "single" : std::to_string(timeWindow).c_str()) ); 

  ///=========== Find all BIN file with runNum 
  std::string pattern = "DataR_run_" + std::to_string(runNum);
  printf("Searching file with pattern %s....\n", pattern.c_str());

  std::vector<std::pair<int, std::string>> inFileName = findDataFiles(pattern);
  int nFile = inFileName.size();

  printf("-------> Out file name : %s \n", outFileName.Data());
  printf("=========================================\n"); 
  printf("       Time Window = %ld ns = %.1f us\n", timeWindow, timeWindow/1000.);
  printf("    Max multiplity = %d hits/event (hard coded)\n", MAX_MULTI);
  printf("========================================= Number of Files : %d \n", nFile);

  for( int i = 0 ; i < nFile; i++ ){
    printf("   %3d | %s\n", i, inFileName[i].second.c_str());
  }

  return 0;

  /*

  unsigned long long int totalHitCount = 0;

  BinReader reader[nFile];
  for(int i = 0; i < nFile; i++ ){
    printf(">>>>>>>> %2d | %s \n", i, inFileName[i].Data());
    reader[i].OpenFile(inFileName[i], 0, true);
    // if( timeOffsetList.size() > 1 ){
    //   for (size_t i = 1 ; timeOffsetList.size(); i++) {
    //     int sn = extractDigiSN(inFileName[i].Data(), timeOffsetList[0].first);
    //     if( sn == timeOffsetList[i].first ) reader[i].SetTimeOffset(timeOffsetList[i].second);
    //   }
    // }
    reader[i].ScanNumHit();
    totalHitCount += reader[i].GetNumHit();
  }

  printf("======================= total Hit Count : %llu\n", totalHitCount);

  printf("========================================= Initializing tree...\n"); 

  unsigned long long                evID = 0;
  unsigned int                     multi = 0;
  unsigned short           sn[MAX_MULTI] = {0}; 
  unsigned short           ch[MAX_MULTI] = {0}; 
  unsigned short            e[MAX_MULTI] = {0};
  // unsigned short           e2[MAX_MULTI] = {0};
  unsigned long long      e_t[MAX_MULTI] = {0};
  // unsigned short          e_f[MAX_MULTI] = {0};
  unsigned short  traceLength[MAX_MULTI];
  short trace[MAX_MULTI][MAX_TRACE_LENGTH];

  TFile * outRootFile = new TFile(outFileName, "recreate");
  TTree * tree = new TTree("tree", outFileName);
  tree->Branch("evID",           &evID, "event_ID/l"); 
  tree->Branch("multi",         &multi, "multi/i"); 
  tree->Branch("sn",                sn, "sn[multi]/s");
  tree->Branch("ch",                ch, "ch[multi]/s");
  tree->Branch("e",                  e, "e[multi]/s");
  // tree->Branch("e2",                e2, "e2[multi]/s");
  tree->Branch("e_t",              e_t, "e_t[multi]/l");
  // tree->Branch("e_f",              e_f, "e_f[multi]/s");
  tree->Branch("traceLength", traceLength, "traceLength[multi]/s");
  
  if( traceOn ) {
    tree->Branch("trace", trace,"trace[multi][MAX_TRACE_LENGTH]/S");
    tree->GetBranch("trace")->SetCompressionSettings(205);
  }

  printf("========================================= Event Building ...\n"); 

  //set TimeWindow to ps;
  timeWindow *= 1000;

  // record the first time stamp and last time stamp
  unsigned long long tStart = 0;
  unsigned long long tEnd = 0;

  unsigned long long t0 = -1;
  short g0 = 0 ;
  int nFileFinished = 0;
  multi = 0;
  evID = 0;

  std::vector<Data> events;
  long long ID[nFile]; // filled hit ID
  for( int i = 0 ; i < nFile ; i++ ) {
    reader[i].ReadBlock();
    SetTimeOffset(&reader[i], timeOffsetList);
    ID[i] = -1;
  }

  unsigned long long hitProcessed = 0;

  do{

    if( debug ) printf("################################ ev build %llu \n", evID);
    events.clear();
    
    // //find earliest time group;
    t0 = -1;
    for( short i = 0; i < nFile; i++){
      if( debug ) printf("NNNNNNNNNNNNNNN %d | ID: %lld, %lld, %lld \n", i, ID[i], reader[i].GetHitID(), reader[i].GetNumHit() );
      if( ID[i] + 1 >= reader[i].GetNumHit())  continue;
      // if( ID[i] >= reader[i].GetHitID() ) reader[i].ReadBlock();
      if( debug ) reader[i].data.Print();
      if( reader[i].data.TimeStamp < t0 ) {
        t0 = reader[i].data.TimeStamp;
        g0 = i;
      }
    }

    if( evID == 0 ) tStart = t0;
    if( debug ) printf("First timestamp is %llu, file ID : %u\n", t0, g0);

    for( short i = 0; i < nFile; i++){

      if( ID[i] + 1 >= reader[i].GetNumHit() ) continue;
      
      if( timeWindow >= 0 ){
        int ret = 0;
        do{
          if( (long int)( reader[i].data.TimeStamp - t0) <= timeWindow ){
            events.push_back(reader[i].data);
            ID[i] ++;
          }else{
            break;
          }
          ret = reader[i].ReadBlock();
          SetTimeOffset(&reader[i], timeOffsetList);
        }while( ret == 1  );
      }else{
        events.push_back(reader[g0].data);
        ID[g0] ++;
        reader[g0].ReadBlock();
        SetTimeOffset(&reader[g0], timeOffsetList);
        break;
      }

      //if( timeWindow < 0) break;
    }

    if( events.size() > 1 ){
      std::sort(events.begin(), events.end(), [](const Data& a, const Data& b) {
        return a.TimeStamp < b.TimeStamp;
      });
    }
    
    tEnd = events.back().TimeStamp;

    hitProcessed += events.size();    
    if( hitProcessed % (traceOn ? 10000 : 10000) == 0 ) printf("hit Porcessed %llu/%llu hit....%.2f%%\n\033[A\r", hitProcessed, totalHitCount, hitProcessed*100./totalHitCount);

    multi = events.size() ;
    if( events.size() >= MAX_MULTI ) {
      printf("\033[31m event %lld has size = %d > MAX_MULTI = %d \033[0m\n", evID, multi, MAX_MULTI);
      // for( int k = 0 ; k < multi; k++){
      //   printf("%d | %d  %2d %llu | %llu , %llu , %llu\n", k, events[k].BoardID, events[k].Channel, events[k].TimeStamp, (long int)( events[k].TimeStamp - t0), t0, timeWindow);
      // }
      multi = MAX_MULTI;
    }
    if( debug ) printf("=================================== filling data | %u \n", multi);

    for( size_t p = 0; p < multi ; p ++ ) {
      if( debug ) {printf("%4zu | ", p); events[p].Print();}

      sn[p] = events[p].BoardID;
      ch[p] = events[p].Channel;
      e[p] = events[p].Energy;
      // e2[p] = events[p].Energy_short;
      e_t[p] = events[p].TimeStamp / 1000;
      // e_f[p] = events[p].TimeStamp % 1000;
      
      traceLength[p] = events[p].NSample;
      if( traceOn ){
        if( traceLength[p] > MAX_TRACE_LENGTH ) {
          printf("\033[31m event %lld has trace length = %d > MAX_TRACE_LENGTH = %d \033[0m\n", evID, traceLength[p], MAX_TRACE_LENGTH);
          traceLength[p] = MAX_TRACE_LENGTH;
        }

        for( int hh = 0; hh < traceLength[p]; hh++){
          trace[p][hh] = events[p].Trace[hh];
        }
      }
    }

    outRootFile->cd();
    tree->Fill();
    // tree->Write();

    multi = 0;
    evID ++;

    nFileFinished = 0;
    for( int i = 0; i < nFile; i++ ){
      // printf(" %d | %lld  %lld \n", i, ID[i], reader[i].GetNumHit() );
      if( ID[i] + 1 >= reader[i].GetNumHit() ) nFileFinished ++ ; 
    }

    if( debug ) printf(" >>>> nFileFinished : %d \n", nFileFinished);

    if( events.size() == 0 ) break;

    if( debug > 1 && hitProcessed > debug ) break;

  }while( nFileFinished < nFile );

  tree->Write();

  unsigned int runEndTime = getTime_us();
  double runTime = (runEndTime - runStartTime) * 1e-6;
  printf("######################################### finished.\n");
  printf(" event building time = %.2f sec = %.2f min\n", runTime, runTime/60.);
  // printf("  total events built = %llu by event builder (%llu in tree)\n", evID, tree->GetEntriesFast());
  printf(" total hit processed = %llu\n", hitProcessed);
  printf("  total events built = %llu by event builder\n", evID);
  double tDuration_sec = (tEnd - tStart) * 1e-12;
  printf("     first timestamp = %20llu ps\n", tStart);
  printf("      last timestamp = %20llu ps\n", tEnd);
  printf(" total data duration = %.2f sec = %.2f min\n", tDuration_sec, tDuration_sec/60.);
  printf("========================================> saved to %s \n", outFileName.Data());

  TMacro info;
  info.AddLine(Form("tStart= %20llu ns",tStart));
  info.AddLine(Form("  tEnd= %20llu ns",tEnd));
  info.Write("info");
  outRootFile->Close();

  return 0;
  */
}