modified: .gitignore

new file:   Armory/EventBuilder.cpp
	modified:   MakeVertex.C
	modified:   run_sx3.sh

.gitignore

@@ -1,4 +1,3 @@
-EventBuilder*
 *.d
 *.so
 *.pcm
@@ -32,3 +31,5 @@ QQQStage2.C
 anasen_fem/scalars.dat.names
 anasen_fem/He96_CO2_4_260Torr.gas
 anasen_fem/heavy_ion_track.csv
+Armory/EventBuilder
+EventBuilder

Armory/EventBuilder.cpp

@@ -0,0 +1,355 @@
+#include "fsuReader.h"
+
+#include "TROOT.h"
+#include "TSystem.h"
+#include "TClonesArray.h"
+#include "TGraph.h"
+#include "TFile.h"
+#include "TTree.h"
+#include "TMacro.h"
+
+
+#define MAX_MULTI  2000
+
+struct FileInfo{
+
+  std::string fileName;
+  int fileID;
+  unsigned long hitCount;
+
+};
+
+//^#############################################################
+//^#############################################################
+int main(int argc, char **argv) {
+  
+  printf("=========================================\n");
+  printf("===      *.fsu Events Builder         ===\n");
+  printf("=========================================\n");  
+  if (argc < 6)    {
+    printf("Incorrect number of arguments:\n");
+    printf("%s [timeWindow] [withTrace] [verbose] [batchSize] [inFile1]  [inFile2] .... \n", argv[0]);
+    printf("    timeWindow : in ns, -1 = no event building \n");   
+    printf("     withTrace : 0 for no trace, 1 for trace \n");   
+    printf("       verbose : > 0 for debug  \n");   
+    printf("     batchSize : the size of hit in a batch \n");   
+    printf("    Output file name is contructed from inFile1 \n");   
+    printf("\n");
+    printf(" Example: %s 0 0 0 10000 '\\ls -1 *001*.fsu'\n", argv[0]);
+    printf("\n\n");
+
+    return 1;
+  }
+
+  uInt runStartTime = getTime_us();
+
+  ///============= read input
+  long timeWindow = atoi(argv[1]);
+  bool traceOn = atoi(argv[2]);
+  unsigned int debug = atoi(argv[3]);
+  unsigned int batchSize = atoi(argv[4]);
+  int nFile = argc - 5;
+  TString inFileName[nFile];
+  for( int i = 0 ; i < nFile ; i++){ inFileName[i] = argv[i+5];}
+  
+  /// Form outFileName;
+  TString outFileName = inFileName[0];
+  int pos = outFileName.Last('/');
+  pos = outFileName.Index("_", pos+1); // find next "_"
+  pos = outFileName.Index("_", pos+1); // find next "_"
+  if( nFile == 1 ) pos = outFileName.Index("_", pos+1); // find next "_", S/N
+  outFileName.Remove(pos); // remove the rest
+  outFileName += "_" + std::to_string(timeWindow);
+  outFileName += ".root";
+  printf("-------> Out file name : %s \n", outFileName.Data());
+  printf("========================================= Number of Files : %d \n", nFile);
+  for( int i = 0; i < nFile; i++) printf("%2d | %s \n", i, inFileName[i].Data());
+  printf("=========================================\n"); 
+  printf("    Time Window = %ld ns = %.1f us\n", timeWindow, timeWindow/1000.);
+  printf("  Include Trace = %s\n", traceOn ? "Yes" : "No");
+  printf("    Debug level = %d\n", debug);
+  printf("     Batch size = %d events/file\n", batchSize);
+  printf(" Max multiplity = %d hits/event (hard coded)\n", MAX_MULTI);
+  printf("========================================= Grouping files\n");  
+
+  std::vector<std::vector<FileInfo>> fileGroupList; // fileName and ID = SN * 1000 + index
+  std::vector<FileInfo> fileList;
+
+  unsigned long long int totalHitCount = 0;
+
+  FSUReader * readerA = new FSUReader(inFileName[0].Data(), 1, 1);
+  readerA->ScanNumBlock(0,0);
+  FileInfo fileInfo = {inFileName[0].Data(), readerA->GetSN() * 1000 +  readerA->GetFileOrder(), readerA->GetHitCount()};
+  fileList.push_back(fileInfo);
+  totalHitCount += readerA->GetHitCount();
+
+  for( int i = 1; i < nFile; i++){
+    FSUReader * readerB = new FSUReader(inFileName[i].Data(), 1, 0);
+    readerB->ScanNumBlock(0,0);
+    totalHitCount += readerB->GetHitCount();
+    fileInfo = {inFileName[i].Data(), readerB->GetSN() * 1000 +  readerB->GetFileOrder(), readerB->GetHitCount()};
+    
+    if( readerA->GetSN() == readerB->GetSN() ){
+      fileList.push_back(fileInfo);
+    }else{
+      fileGroupList.push_back(fileList);
+      fileList.clear();
+      fileList.push_back(fileInfo);
+    }
+  
+    delete readerA;
+    readerA = readerB;
+  }
+  fileGroupList.push_back(fileList);
+  delete readerA;
+
+  printf("======================= total Hit Count : %llu\n", totalHitCount);
+
+  for( size_t i = 0; i < fileGroupList.size(); i++){
+    printf("group ----- %ld \n", i);
+
+    //sort by ID
+    std::sort(fileGroupList[i].begin(), fileGroupList[i].end(), [](const FileInfo & a, const FileInfo & b) {
+      return a.fileID < b.fileID;
+    });
+
+    for( size_t j = 0; j < fileGroupList[i].size(); j++){
+      printf("%3ld | %8d | %9lu| %s \n", j, fileGroupList[i][j].fileID, fileGroupList[i][j].hitCount, fileGroupList[i][j].fileName.c_str() );
+    }
+
+  }
+
+  // //*====================================== create tree
+  TFile * outRootFile = new TFile(outFileName, "recreate");
+  TTree * tree = new TTree("tree", outFileName);
+
+  unsigned long long                evID = 0;
+  unsigned int                     multi = 0;
+  unsigned short           sn[MAX_MULTI] = {0}; /// board SN
+  unsigned short           ch[MAX_MULTI] = {0}; /// chID
+  unsigned short            e[MAX_MULTI] = {0}; /// 15 bit
+  unsigned short           e2[MAX_MULTI] = {0}; /// 15 bit
+  unsigned long long      e_t[MAX_MULTI] = {0}; /// timestamp 47 bit
+  unsigned short          e_f[MAX_MULTI] = {0}; /// fine time 10 bit 
+  unsigned short  traceLength[MAX_MULTI];
+
+  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_timestamp[multi]/l");
+  tree->Branch("e_f",              e_f, "e_fineTime[multi]/s");
+  tree->Branch("traceLength", traceLength, "traceLength[multi]/s");
+  
+  TClonesArray * arrayTrace = nullptr;
+  TGraph * trace = nullptr;
+
+  if( traceOn ) {
+    arrayTrace = new TClonesArray("TGraph");
+    tree->Branch("trace",       arrayTrace, 2560000);
+    arrayTrace->BypassStreamer();
+  }
+
+  //*======================================= Open files
+  printf("========================================= Open files & Build Events.\n"); 
+
+  const short nGroup = fileGroupList.size();
+  std::vector<Hit> hitList[nGroup];
+  
+  FSUReader ** reader = new FSUReader * [nGroup];
+  ulong ID[nGroup];
+  for( short i = 0; i < nGroup; i++){
+    std::vector<std::string> fList;
+    for( size_t j = 0; j < fileGroupList[i].size(); j++){
+      fList.push_back( fileGroupList[i][j].fileName );
+    }
+    reader[i] = new FSUReader(fList, 600, debug);
+    hitList[i] = reader[i]->ReadBatch(batchSize, debug );
+    reader[i]->PrintHitListInfo(&hitList[i], "hitList-" + std::to_string(reader[i]->GetSN()));
+    ID[i] = 0;
+    if( debug ) {
+
+      for( size_t p = 0; p < 10; p ++ ){
+        if( hitList[i].size() <= p ) break;
+        hitList[i][p].Print();
+      }
+
+    }
+  }
+
+  unsigned long long tStart = 0;
+  unsigned long long tEnd = 0;
+
+  //find earliest time group;
+  unsigned long long t0 = -1;
+  short g0 = 0 ;
+  for( short i = 0; i < nGroup; i++){
+    if( hitList[i].size() == 0 ) continue;
+    if( hitList[i][0].timestamp < t0 ) {
+      t0 = hitList[i][0].timestamp;
+      g0 = i;
+    }
+  }
+  tStart = t0;
+  if( debug ) printf("First timestamp is %llu, group : %u\n", t0, g0);
+
+  int nFileFinished = 0;
+  multi = 0;
+  evID = 0;
+  std::vector<Hit> events;
+
+  unsigned long long hitProcessed = 0;
+
+  do{
+
+    //*============= Build events from hitList[i]
+    if( debug ) printf("################################ ev build %llu \n", evID);
+    events.clear();
+
+    for( short i = 0; i < nGroup; i++){
+      short ig = (i + g0 ) % nGroup;
+
+      if( hitList[ig].size() == 0 ) continue;
+
+      //chekc if reached the end of hitList
+      if( ID[ig] >= hitList[ig].size() ) {
+        hitList[ig] = reader[ig]->ReadBatch(batchSize, debug + 1);
+        if( debug ) reader[ig]->PrintHitListInfo( &hitList[ig], "hitList-" + std::to_string(ig));
+        ID[ig] = 0;
+        if( hitList[ig].size() == 0 ) continue;
+      }
+
+      if( timeWindow >= 0 ){
+      
+        do{
+
+          if( (long int)(hitList[ig].at(ID[ig]).timestamp  - t0) <= timeWindow ){
+            events.push_back(hitList[ig].at(ID[ig]));
+            ID[ig] ++;
+          }else{
+            break;
+          }
+
+          //check if reached the end of hitList
+          if( ID[ig] >= hitList[ig].size() ) {
+            hitList[ig] = reader[ig]->ReadBatch(batchSize, debug);
+            if( debug ) reader[ig]->PrintHitListInfo( &hitList[ig], "hitList-" + std::to_string(ig));
+            ID[ig] = 0;
+            if( hitList[ig].size() == 0 ) break;
+          }
+      
+        }while(ID[ig] < hitList[ig].size());
+
+      }else{
+        events.push_back(hitList[ig].at(ID[ig]));
+        ID[ig] ++;
+      }
+
+      if( timeWindow < 0) break;
+
+    }
+
+    if( events.size() > 1 ){
+      std::sort(events.begin(), events.end(), [](const Hit& a, const Hit& b) {
+        return a.timestamp < b.timestamp;
+      });
+    }
+    
+    tEnd = events.back().timestamp;
+
+    hitProcessed += events.size();    
+    if( hitProcessed % 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("event %lld has size = %d > MAX_MULTI = %d\n", evID, multi, MAX_MULTI);
+      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].sn;
+      ch[p] = events[p].ch;
+      e[p] = events[p].energy;
+      e2[p] = events[p].energy2;
+      e_t[p] = events[p].timestamp;
+      e_f[p] = events[p].fineTime;
+      
+      traceLength[p] = events[p].traceLength;
+      if( traceOn ){
+        trace = (TGraph *) arrayTrace->ConstructedAt(multi, "C");
+        trace->Clear();
+        for( int hh = 0; hh < traceLength[multi]; hh++){
+          trace->SetPoint(hh, hh, events[p].trace[hh]);
+        }
+      }
+    }
+
+    outRootFile->cd();
+    tree->Fill();
+    // tree->Write();
+    multi = 0;
+    evID ++;
+
+    ///===================== find the next first timestamp
+    t0 = -1;
+    g0 = -1;
+
+    for( int i = 0; i < nGroup; i++) {
+      if( hitList[i].size() == 0 ) continue;
+      if( hitList[i][ID[i]].timestamp < t0 ) {
+        t0 = hitList[i][ID[i]].timestamp;
+        g0 = i;
+      }
+    } 
+    if( debug ) printf("Next First timestamp is %llu, group : %u\n", t0, g0);
+
+    //*=============
+    nFileFinished = 0;
+    for( int i = 0 ; i < nGroup; i++) {
+      if( hitList[i].size() == 0 ) {
+        nFileFinished ++;
+        continue;
+      }else{
+        if( ID[i] >= hitList[i].size( )) {
+          hitList[i] = reader[i]->ReadBatch(batchSize, debug);
+          ID[i] = 0;
+          if( hitList[i].size() == 0 ) nFileFinished ++;
+        }
+      }
+    }
+    if( debug > 1 ) printf("========== nFileFinished : %d\n", nFileFinished);
+
+  }while( nFileFinished < nGroup);
+
+  tree->Write();
+
+  uInt 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());
+  double tDuration_sec = (tEnd - tStart) * 1e-9;
+  printf("     first timestamp = %20llu ns\n", tStart);
+  printf("      last timestamp = %20llu ns\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();
+  
+  for( int i = 0; i < nGroup; i++) delete reader[i];
+  delete [] reader;
+
+  return 0;
+
+}
+

MakeVertex.C

@@ -1145,7 +1145,8 @@ Bool_t MakeVertex::Process(Long64_t entry)
             plotter->Fill2D("pcz_vs_sx3z_strip="+std::to_string(sx3event.ch2)+"_a"+std::to_string(pcevent.multi1)+"_c"+std::to_string(pcevent.multi2),300,0,100,600,-200,200,sx3z,pcevent.pos.Z(),"sx3_vs_pc_zcorr");
 
             plotter->Fill2D("pcdEC_vs_sx3z_strip="+std::to_string(sx3event.ch2)+"_a"+std::to_string(pcevent.multi1)+"_c"+std::to_string(pcevent.multi2),800,0,20000,600,-200,200,pcevent.Energy2,sx3z,"sx3_vs_pc_zcorr");
-            plotter->Fill2D("pcdEA_vs_sx3z_strip="+std::to_string(sx3event.ch2)+"_a"+std::to_string(pcevent.multi1)+"_c"+std::to_string(pcevent.multi2),800,0,20000,600,-200,200,pcevent.Energy1,sx3z,"sx3_vs_pc_zcorr");*/
+            plotter->Fill2D("pcdEA_vs_sx3z_strip="+std::to_string(sx3event.ch2)+"_a"+std::to_string(pcevent.multi1)+"_c"+std::to_string(pcevent.multi2),800,0,20000,600,-200,200,pcevent.Energy1,sx3z,"sx3_vs_pc_zcorr");
+            */
 
             /*for(auto cc: cClusters)
              for(auto ac: aClusters) {

run_sx3.sh

@@ -1,5 +1,5 @@
 #Alpha runs at different spacer positions
-rm results_run*.root
+# rm results_run*.root
 export flipa=0
 export anode_offset=0
 export DATASET="27Al"
@@ -22,14 +22,14 @@ export flip180="0"
 #root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_009_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run09.root;
 if [[ 1 -eq 1 ]]; then
 #export timecut_low=230.0;
-# export timecut_low=400.0;
+export timecut_low=400.0;
 #export timecut_high=400.0;
 #unset timecut_low, timecut_high
 #export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_009_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run09.root;
 #export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_010_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run10.root;
 #export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_011_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run11.root;
 export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_012_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run12.root;
-# exit
+exit
 #export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_013_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run13.root;
 #exit
 fi