rtang/SOLARIS_Analysis
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

update Event Builder

Browse Source
This commit is contained in:
Ryan Tang 2023-09-26 15:45:46 -04:00
parent c9d661e709
commit 8ce25469d9
9 changed files with 479 additions and 206 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
.gitignore vendored
View File

@ -6,6 +6,8 @@ EventBuilder
*.d *.d
*.so *.so
.gdb_history
data_raw data_raw
root_data root_data

4
.vscode/settings.json vendored
View File

@ -119,7 +119,9 @@
"locale": "cpp", "locale": "cpp",
"__verbose_abort": "cpp", "__verbose_abort": "cpp",
"Monitors_Util.C": "cpp", "Monitors_Util.C": "cpp",
"Monitor_Util.C": "cpp" "Monitor_Util.C": "cpp",
"script_single.C": "cpp",
"script_multi.C": "cpp"
}, },
"better-comments.multilineComments": true, "better-comments.multilineComments": true,

37
README.md Normal file
View File

@ -0,0 +1,37 @@
# Introduction
This is the analysis package for the SOLARIS DAQ. It is supposed to be the analysis path for the SOLARIS DAQ.
The folder struture is
Analysis
├── README.md
├── SetupNewExp // bash script to create new branch and raw data folder
├── SOLARIS.sh // bash script to define some env variable and functions
├── armory // analysis codes, independent from experiment.
├── Cleopatra // Swaper for DWBA code Ptolomey
├── data_raw // should be the symbolic link to the raw data, created by SetUpNewExp
├── root_data // symbolic link to converted root file, created by SetUpNewExp
└── working // working directory, depends on experiment.
# Event Builder
The EventBuilder is at the armory. It depends on the Hit.h and SolReader.h.
## Hit.h
The Hit class stores a hit (or a data block)
## SolReader.h
The SolReader class read the sol file. It can be loaded in CERN ROOT alone.

104
armory/EventBuilder.cpp
View File

@ -10,11 +10,11 @@
//#include "TClonesArray.h" // plan to save trace as TVector with TClonesArray //#include "TClonesArray.h" // plan to save trace as TVector with TClonesArray
//#include "TVector.h" //#include "TVector.h"
#define MAX_ID 64 #define MAX_MULTI 64
#define MAX_TRACE_LEN 2500 #define MAX_TRACE_LEN 2500
SolReader ** reader; SolReader ** reader;
Event ** evt; Hit ** hit;
unsigned long totFileSize = 0; unsigned long totFileSize = 0;
unsigned long processedFileSize = 0; unsigned long processedFileSize = 0;
@ -25,48 +25,51 @@ std::vector<std::vector<int>> group; // group[i][j], i = group ID, j = group mem
void findEarliestTime(int &fileID, int & digiID){ void findEarliestTime(int &fileID, int & digiID){
unsigned long firstTime = 0; unsigned long firstTime = 0;
for( int i = 0; i < (int) activeFileID.size(); i++){ for( int i = 0; i < (int) activeFileID.size(); i++){
int id = activeFileID[i]; int id = activeFileID[i];
if( i == 0 ) { if( i == 0 ) {
firstTime = evt[id]->timestamp; firstTime = hit[id]->timestamp;
fileID = id; fileID = id;
digiID = i; //printf("%d | %d %lu %d | %d \n", id, reader[id]->GetBlockID(), hit[id]->timestamp, hit[id]->channel, (int) activeFileID.size());
//printf("%d | %ld %lu %d | %d \n", id, reader[id]->GetBlockID(), evt[id]->timestamp, evt[id]->channel, (int) activeFileID.size());
continue; continue;
} }
if( evt[id]->timestamp <= firstTime) { if( hit[id]->timestamp <= firstTime) {
firstTime = evt[id]->timestamp; firstTime = hit[id]->timestamp;
fileID = id; fileID = id;
digiID = i; //printf("%d | %d %lu %d | %d \n", id, reader[id]->GetBlockID(), hit[id]->timestamp, hit[id]->channel, (int) activeFileID.size());
//printf("%d | %ld %lu %d | %d \n", id, reader[id]->GetBlockID(), evt[id]->timestamp, evt[id]->channel, (int) activeFileID.size());
} }
} }
} }
unsigned long long evID = 0; unsigned long long evID = 0;
int multi = 0; unsigned int multi = 0;
int bd[MAX_ID] = {0}; unsigned short bd[MAX_MULTI] = {0};
int ch[MAX_ID] = {0}; unsigned short sn[MAX_MULTI] = {0};
int e[MAX_ID] = {0}; unsigned short ch[MAX_MULTI] = {0};
unsigned long long e_t[MAX_ID] = {0}; unsigned short e[MAX_MULTI] = {0};
unsigned short lowFlag[MAX_ID] = {0}; unsigned short e2[MAX_MULTI] = {0}; //for PSD energy short
unsigned short highFlag[MAX_ID] = {0}; unsigned long long e_t[MAX_MULTI] = {0};
int traceLen[MAX_ID] = {0}; unsigned short lowFlag[MAX_MULTI] = {0};
int trace[MAX_ID][MAX_TRACE_LEN] = {0}; unsigned short highFlag[MAX_MULTI] = {0};
int traceLen[MAX_MULTI] = {0};
int trace[MAX_MULTI][MAX_TRACE_LEN] = {0};
void fillData(int &fileID, int &digiID, const bool &saveTrace){ void fillData(int &fileID, int &digiID, const bool &saveTrace){
bd[multi] = digiID; bd[multi] = digiID;
ch[multi] = evt[fileID]->channel; sn[multi] = digiID;
e[multi] = evt[fileID]->energy; ch[multi] = hit[fileID]->channel;
e_t[multi] = evt[fileID]->timestamp; e[multi] = hit[fileID]->energy;
lowFlag[multi] = evt[fileID]->flags_low_priority; e2[multi] = hit[fileID]->energy_short;
highFlag[multi] = evt[fileID]->flags_high_priority; e_t[multi] = hit[fileID]->timestamp;
lowFlag[multi] = hit[fileID]->flags_low_priority;
highFlag[multi] = hit[fileID]->flags_high_priority;
if( saveTrace ){ if( saveTrace ){
traceLen[multi] = evt[fileID]->traceLenght; traceLen[multi] = hit[fileID]->traceLenght;
for( int i = 0; i < TMath::Min(traceLen[multi], MAX_TRACE_LEN); i++){ for( int i = 0; i < TMath::Min(traceLen[multi], MAX_TRACE_LEN); i++){
trace[multi][i] = evt[fileID]->analog_probes[0][i]; trace[multi][i] = hit[fileID]->analog_probes[0][i];
} }
} }
@ -114,14 +117,11 @@ int main(int argc, char ** argv){
//*======================================== setup reader //*======================================== setup reader
reader = new SolReader*[nFile]; reader = new SolReader*[nFile];
evt = new Event *[nFile]; hit = new Hit *[nFile];
for( int i = 0 ; i < nFile ; i++){ for( int i = 0 ; i < nFile ; i++){
reader[i] = new SolReader(inFileName[i].Data()); reader[i] = new SolReader(inFileName[i].Data());
evt[i] = reader[i]->evt; //TODO check is file open propertly hit[i] = reader[i]->hit; //TODO check is file open propertly
//reader[i]->ScanNumBlock();
reader[i]->ReadNextBlock(); // read the first block reader[i]->ReadNextBlock(); // read the first block
} }
@ -144,16 +144,20 @@ int main(int argc, char ** argv){
int digiID = f1.Remove(pos).Atoi(); int digiID = f1.Remove(pos).Atoi();
fn.Remove(0, pos+1); fn.Remove(0, pos+1);
pos = fn.Last('_'); // remove digi serial num pos = fn.Last('_'); // digi serial num
f1 = fn;
int digisn = f1.Remove(pos).Atoi();
fn.Remove(0, pos+1); fn.Remove(0, pos+1);
pos = fn.First('.'); // get the file id; pos = fn.First('.'); // get the file id;
int fileID = fn.Remove(pos).Atoi(); int indexID = fn.Remove(pos).Atoi();
std::vector<int> haha = {i, digiID, fileID}; int fileID = i;
std::vector<int> haha = {fileID, digiID, indexID, digisn};
idList.push_back(haha); idList.push_back(haha);
} }
// sort by digiID
std::sort(idList.begin(), idList.end(), [](const std::vector<int>& a, const std::vector<int>& b){ std::sort(idList.begin(), idList.end(), [](const std::vector<int>& a, const std::vector<int>& b){
if (a[1] == b[1]) { if (a[1] == b[1]) {
return a[2] < b[2]; return a[2] < b[2];
@ -161,11 +165,11 @@ int main(int argc, char ** argv){
return a[1] < b[1]; return a[1] < b[1];
}); });
group.clear(); group.clear(); // group[i][j], i is the group Index = digiID
int last_id = 0; int last_id = 0;
std::vector<int> kaka; std::vector<int> kaka;
for( int i = 0; i < (int) idList.size() ; i++){ for( int i = 0; i < (int) idList.size() ; i++){
if( i == 0 ) { if( i == 0 ) {
kaka.clear(); kaka.clear();
last_id = idList[i][1]; last_id = idList[i][1];
kaka.push_back(idList[i][0]); kaka.push_back(idList[i][0]);
@ -208,10 +212,12 @@ int main(int argc, char ** argv){
TTree * tree = new TTree("tree", outFileName); TTree * tree = new TTree("tree", outFileName);
tree->Branch("evID", &evID, "event_ID/l"); tree->Branch("evID", &evID, "event_ID/l");
tree->Branch("multi", &multi, "multi/I"); tree->Branch("multi", &multi, "multi/i");
tree->Branch("bd", bd, "board[multi]/I"); tree->Branch("bd", bd, "board[multi]/s");
tree->Branch("ch", ch, "channel[multi]/I"); tree->Branch("sn", sn, "sn[multi]/s");
tree->Branch("e", e, "energy[multi]/I"); tree->Branch("ch", ch, "channel[multi]/s");
tree->Branch("e", e, "energy[multi]/s");
tree->Branch("e2", e2, "energy_short[multi]/s");
tree->Branch("e_t", e_t, "energy_timestamp[multi]/l"); tree->Branch("e_t", e_t, "energy_timestamp[multi]/l");
tree->Branch("lowFlag", lowFlag, "lowFlag[multi]/s"); tree->Branch("lowFlag", lowFlag, "lowFlag[multi]/s");
tree->Branch("highFlag", highFlag, "highFlag[multi]/s"); tree->Branch("highFlag", highFlag, "highFlag[multi]/s");
@ -241,15 +247,16 @@ int main(int argc, char ** argv){
findEarliestTime(fileID, digiID); findEarliestTime(fileID, digiID);
fillData(fileID, digiID, saveTrace); fillData(fileID, digiID, saveTrace);
unsigned long firstTimeStamp = evt[fileID]->timestamp; unsigned long firstTimeStamp = hit[fileID]->timestamp;
unsigned long lastTimeStamp = 0; unsigned long lastTimeStamp = 0;
int last_precentage = 0; int last_precentage = 0;
while((activeFileID.size() > 0)){ while((activeFileID.size() > 0)){
findEarliestTime(fileID, digiID); findEarliestTime(fileID, digiID);
if( evt[fileID]->timestamp - e_t[0] < timeWindow ){ if( reader[fileID]->IsEndOfFile() ) relay[digiID] = true;
if( hit[fileID]->timestamp - e_t[0] < timeWindow ){
fillData(fileID, digiID, saveTrace); fillData(fileID, digiID, saveTrace);
}else{ }else{
outRootFile->cd(); outRootFile->cd();
@ -259,7 +266,7 @@ int main(int argc, char ** argv){
multi = 0; multi = 0;
fillData(fileID, digiID, saveTrace); fillData(fileID, digiID, saveTrace);
} }
///========= check is file finished. ///========= check is file finished.
if( relay[digiID]){ if( relay[digiID]){
groupIndex[digiID] ++; groupIndex[digiID] ++;
@ -271,8 +278,7 @@ int main(int argc, char ** argv){
} }
relay[digiID] = false; relay[digiID] = false;
} }
if( reader[fileID]->IsEndOfFile() ) relay[digiID] = true;
///========= calculate progress ///========= calculate progress
processedFileSize = 0; processedFileSize = 0;
for( int p = 0; p < (int) group.size(); p ++){ for( int p = 0; p < (int) group.size(); p ++){
@ -302,7 +308,7 @@ int main(int argc, char ** argv){
double percentage = processedFileSize * 100/ totFileSize; double percentage = processedFileSize * 100/ totFileSize;
printf("Processed : %llu, %.0f%% | %lu/%lu \n", evID, percentage, processedFileSize, totFileSize); printf("Processed : %llu, %.0f%% | %lu/%lu \n", evID, percentage, processedFileSize, totFileSize);
lastTimeStamp = evt[fileID]->timestamp; lastTimeStamp = hit[fileID]->timestamp;
//*=========================================== save file //*=========================================== save file
outRootFile->cd(); outRootFile->cd();
tree->Fill(); tree->Fill();

148
armory/Event.h → armory/Hit.h
View File

@ -1,5 +1,5 @@
#ifndef EVENT_H #ifndef HIT_H
#define EVENT_H #define HIT_H
#include <stdio.h> #include <stdio.h>
#include <cstdlib> #include <cstdlib>
@ -8,27 +8,46 @@
#define MaxTraceLenght 8100 #define MaxTraceLenght 8100
class Event { enum DataFormat{
ALL = 0,
OneTrace = 1,
NoTrace = 2,
Minimum = 3,
RAW = 0x0A,
};
namespace DPPType{
const std::string PHA = "DPP_PHA";
const std::string PSD = "DPP_PSD";
};
class Hit {
public: public:
unsigned short dataType; unsigned short dataType;
std::string DPPType;
///============= for dpp-pha ///============= for dpp-pha
uint8_t channel; // 6 bit uint8_t channel; // 6 bit
uint16_t energy; // 16 bit uint16_t energy; // 16 bit
uint64_t timestamp; // 48 bit uint16_t energy_short; // 16 bit, only for PSD
uint64_t timestamp; // 48 bit
uint16_t fine_timestamp; // 16 bit uint16_t fine_timestamp; // 16 bit
uint16_t flags_low_priority; // 12 bit uint16_t flags_low_priority; // 12 bit
uint16_t flags_high_priority; // 8 bit uint16_t flags_high_priority; // 8 bit
size_t traceLenght; // 64 bit size_t traceLenght; // 64 bit
uint8_t downSampling; // 8 bit uint8_t downSampling; // 8 bit
bool board_fail; bool board_fail;
bool flush; bool flush;
uint8_t analog_probes_type[2]; // 3 bit uint8_t analog_probes_type[2]; // 3 bit for PHA, 4 bit for PSD
uint8_t digital_probes_type[4]; // 4 bit uint8_t digital_probes_type[4]; // 4 bit for PHA, 5 bit for PSD
int32_t * analog_probes[2]; // 18 bit int32_t * analog_probes[2]; // 18 bit
uint8_t * digital_probes[4]; // 1 bit uint8_t * digital_probes[4]; // 1 bit
uint16_t trigger_threashold; // 16 bit uint16_t trigger_threashold; // 16 bit
size_t event_size; // 64 bit size_t event_size; // 64 bit
uint32_t aggCounter; // 32 bit uint32_t aggCounter; // 32 bit
@ -39,17 +58,21 @@ class Event {
bool isTraceAllZero; bool isTraceAllZero;
Event(){ Hit(){
Init(); Init();
} }
~Event(){ ~Hit(){
ClearMemory(); ClearMemory();
} }
void Init(){ void Init(){
DPPType = DPPType::PHA;
dataType = DataFormat::ALL;
channel = 0; channel = 0;
energy = 0; energy = 0;
energy_short = 0;
timestamp = 0; timestamp = 0;
fine_timestamp = 0; fine_timestamp = 0;
downSampling = 0; downSampling = 0;
@ -92,11 +115,12 @@ class Event {
isTraceAllZero = true; isTraceAllZero = true;
} }
void SetDataType(unsigned int type){ void SetDataType(unsigned int type, std::string dppType){
dataType = type; dataType = type;
DPPType = dppType;
ClearMemory(); ClearMemory();
if( dataType == 0xF){ if( dataType == DataFormat::RAW){
data = new uint8_t[20*1024*1024]; data = new uint8_t[20*1024*1024];
}else{ }else{
analog_probes[0] = new int32_t[MaxTraceLenght]; analog_probes[0] = new int32_t[MaxTraceLenght];
@ -128,34 +152,68 @@ class Event {
} }
void PrintEnergyTimeStamp(){ void PrintEnergyTimeStamp(){
printf("ch: %2d, energy: %u, timestamp: %llu ch, traceLenght: %lu\n", channel, energy, timestamp, traceLenght); printf("ch: %2d, energy: %u, timestamp: %lu ch, traceLenght: %lu\n", channel, energy, timestamp, traceLenght);
} }
std::string AnaProbeType(uint8_t probeType){ std::string AnaProbeType(uint8_t probeType){
switch(probeType){
case 0: return "ADC"; if( DPPType == DPPType::PHA){
case 1: return "Time filter"; switch(probeType){
case 2: return "Energy filter"; case 0: return "ADC";
default : return "none"; case 1: return "Time filter";
case 2: return "Energy filter";
default : return "none";
}
}else if (DPPType == DPPType::PSD){
switch(probeType){
case 0: return "ADC";
case 9: return "Baseline";
case 10: return "CFD";
default : return "none";
}
}else{
return "none";
} }
} }
std::string DigiProbeType(uint8_t probeType){ std::string DigiProbeType(uint8_t probeType){
switch(probeType){
case 0: return "Trigger"; if( DPPType == DPPType::PHA){
case 1: return "Time filter armed"; switch(probeType){
case 2: return "Re-trigger guard"; case 0: return "Trigger";
case 3: return "Energy filter baseline freeze"; case 1: return "Time filter armed";
case 4: return "Energy filter peaking"; case 2: return "Re-trigger guard";
case 5: return "Energy filter peaking ready"; case 3: return "Energy filter baseline freeze";
case 6: return "Energy filter pile-up guard"; case 4: return "Energy filter peaking";
case 7: return "Event pile-up"; case 5: return "Energy filter peaking ready";
case 8: return "ADC saturation"; case 6: return "Energy filter pile-up guard";
case 9: return "ADC saturation protection"; case 7: return "Event pile-up";
case 10: return "Post-saturation event"; case 8: return "ADC saturation";
case 11: return "Energy filter saturation"; case 9: return "ADC saturation protection";
case 12: return "Signal inhibit"; case 10: return "Post-saturation event";
default : return "none"; case 11: return "Energy filter saturation";
case 12: return "Signal inhibit";
default : return "none";
}
}else if (DPPType == DPPType::PSD){
switch(probeType){
case 0: return "Trigger";
case 1: return "CFD Filter Armed";
case 2: return "Re-trigger guard";
case 3: return "ADC Input Baseline freeze";
case 20: return "ADC Input OverThreshold";
case 21: return "Charge Ready";
case 22: return "Long Gate";
case 7: return "Pile-Up Trig.";
case 24: return "Short Gate";
case 25: return "Energy Saturation";
case 26: return "Charge over-range";
case 27: return "ADC Input Neg. OverThreshold";
default : return "none";
}
}else{
return "none";
} }
} }
@ -177,9 +235,19 @@ class Event {
//TODO LowPriority //TODO LowPriority
void PrintAll(){ void PrintAll(){
printf("============= Type : %u\n", dataType);
switch(dataType){
case DataFormat::ALL : printf("============= Type : ALL\n"); break;
case DataFormat::OneTrace : printf("============= Type : OneTrace\n"); break;
case DataFormat::NoTrace : printf("============= Type : NoTrace\n"); break;
case DataFormat::Minimum : printf("============= Type : Minimum\n"); break;
case DataFormat::RAW : printf("============= Type : RAW\n"); return; break;
default : return;
}
printf("ch : %2d (0x%02X), fail: %d, flush: %d\n", channel, channel, board_fail, flush); printf("ch : %2d (0x%02X), fail: %d, flush: %d\n", channel, channel, board_fail, flush);
printf("energy: %u, timestamp: %llu, fine_timestamp: %u \n", energy, timestamp, fine_timestamp); if( DPPType == DPPType::PHA ) printf("energy: %u, timestamp: %lu, fine_timestamp: %u \n", energy, timestamp, fine_timestamp);
if( DPPType == DPPType::PSD ) printf("energy: %u, energy_S : %u, timestamp: %lu, fine_timestamp: %u \n", energy, energy_short, timestamp, fine_timestamp);
printf("flag (high): 0x%02X, (low): 0x%03X, traceLength: %lu\n", flags_high_priority, flags_low_priority, traceLenght); printf("flag (high): 0x%02X, (low): 0x%03X, traceLength: %lu\n", flags_high_priority, flags_low_priority, traceLenght);
printf("Agg counter : %u, trigger Thr.: %u, downSampling: %u \n", aggCounter, trigger_threashold, downSampling); printf("Agg counter : %u, trigger Thr.: %u, downSampling: %u \n", aggCounter, trigger_threashold, downSampling);
printf("AnaProbe Type: %s(%u), %s(%u)\n", AnaProbeType(analog_probes_type[0]).c_str(), analog_probes_type[0], printf("AnaProbe Type: %s(%u), %s(%u)\n", AnaProbeType(analog_probes_type[0]).c_str(), analog_probes_type[0],

5
armory/Makefile
View File

@ -1,10 +1,11 @@
CC=g++ CC=g++
CFLAG= -g
ROOTFLAG=`root-config --cflags --glibs` ROOTFLAG=`root-config --cflags --glibs`
all: EventBuilder all: EventBuilder
EventBuilder: EventBuilder.cpp SolReader.h EventBuilder: EventBuilder.cpp SolReader.h Hit.h
$(CC) EventBuilder.cpp -o EventBuilder ${ROOTFLAG} $(CC) $(CFLAG) EventBuilder.cpp -o EventBuilder ${ROOTFLAG}
clean: clean:

167
armory/SolReader.h
View File

@ -8,7 +8,7 @@
#include <unistd.h> #include <unistd.h>
#include <time.h> // time in nano-sec #include <time.h> // time in nano-sec
#include "Event.h" // this is a symblic link to SOLARIS_QT6_DAQ/Event.h #include "Hit.h"
#define tick2ns 8 // 1 tick = 8 ns #define tick2ns 8 // 1 tick = 8 ns
@ -20,7 +20,7 @@ class SolReader {
unsigned int totNumBlock; unsigned int totNumBlock;
unsigned short blockStartIdentifier; unsigned short blockStartIdentifier;
long blockID; unsigned int numBlock;
bool isScanned; bool isScanned;
void init(); void init();
@ -29,7 +29,7 @@ class SolReader {
public: public:
SolReader(); SolReader();
SolReader(std::string fileName, unsigned short dataType = 0); // dataType can auto determine from the data, but remove it will crash.... SolReader(std::string fileName, unsigned short dataType);
~SolReader(); ~SolReader();
void OpenFile(std::string fileName); void OpenFile(std::string fileName);
@ -38,25 +38,25 @@ class SolReader {
void ScanNumBlock(); void ScanNumBlock();
long GetBlockID() const {return blockID;} bool IsEndOfFile() const {return (filePos >= inFileSize ? true : false);}
unsigned int GetBlockID() const {return numBlock - 1;}
unsigned int GetNumBlock() const {return numBlock;}
unsigned int GetTotalNumBlock() const {return totNumBlock;} unsigned int GetTotalNumBlock() const {return totNumBlock;}
unsigned int GetFilePos() const {return filePos;} unsigned int GetFilePos() const {return filePos;}
unsigned int GetFileSize() const {return inFileSize;} unsigned int GetFileSize() const {return inFileSize;}
bool IsEndOfFile() {return (filePos >= inFileSize ? true : false);}
void RewindFile(); void RewindFile();
Event * evt; Hit * hit;
}; };
void SolReader::init(){ void SolReader::init(){
inFileSize = 0; inFileSize = 0;
blockID = -1; numBlock = 0;
filePos = 0; filePos = 0;
totNumBlock = 0; totNumBlock = 0;
evt = new Event(); hit = new Hit();
isScanned = false; isScanned = false;
@ -68,16 +68,15 @@ SolReader::SolReader(){
init(); init();
} }
SolReader::SolReader(std::string fileName, unsigned short dataType){ SolReader::SolReader(std::string fileName, unsigned short dataType = 0){
init(); init();
OpenFile(fileName); OpenFile(fileName);
evt->SetDataType(dataType); hit->SetDataType(dataType, DPPType::PHA);
} }
SolReader::~SolReader(){ SolReader::~SolReader(){
//printf("%s\n", __func__);
if( !inFile ) fclose(inFile); if( !inFile ) fclose(inFile);
delete evt; delete hit;
} }
inline void SolReader::OpenFile(std::string fileName){ inline void SolReader::OpenFile(std::string fileName){
@ -99,135 +98,139 @@ inline int SolReader::ReadBlock(unsigned int index, bool verbose){
if( verbose ) printf("Block index: %u, File Pos: %u byte\n", index, blockPos[index]); if( verbose ) printf("Block index: %u, File Pos: %u byte\n", index, blockPos[index]);
fseek(inFile, blockPos[index], SEEK_CUR); fseek(inFile, blockPos[index], SEEK_CUR);
filePos = blockPos[index]; filePos = blockPos[index];
blockID = index; numBlock = index;
blockID --;
return ReadNextBlock(); return ReadNextBlock();
} }
inline int SolReader::ReadNextBlock(int isSkip){ inline int SolReader::ReadNextBlock(int isSkip){
if( inFile == NULL ) return -1; if( inFile == NULL ) return -1;
if( feof(inFile) ) return -1; if( feof(inFile) ) return -1;
if( filePos >= inFileSize) return -1; if( filePos >= inFileSize) return -1;
fread(&blockStartIdentifier, 2, 1, inFile); fread(&blockStartIdentifier, 2, 1, inFile);
if( (blockStartIdentifier & 0xAAA0) != 0xAAA0 ) { if( (blockStartIdentifier & 0xAA00) != 0xAA00 ) {
printf("header fail.\n"); printf("header fail.\n");
return -2 ; return -2 ;
} }
if( ( blockStartIdentifier & 0xF ) == 15 ){ if( ( blockStartIdentifier & 0xF ) == DataFormat::RAW ){
evt->SetDataType(15); hit->SetDataType(DataFormat::RAW, ((blockStartIdentifier >> 1) & 0xF) == 0 ? DPPType::PHA : DPPType::PSD);
} }
evt->dataType = blockStartIdentifier & 0xF; hit->dataType = blockStartIdentifier & 0xF;
hit->DPPType = ((blockStartIdentifier >> 1) & 0xF) == 0 ? DPPType::PHA : DPPType::PSD;
if( evt->dataType == 0){ //======== same as the dataFormat in Digitizer if( hit->dataType == DataFormat::ALL){
if( isSkip == 0 ){ if( isSkip == 0 ){
fread(&evt->channel, 1, 1, inFile); fread(&hit->channel, 1, 1, inFile);
fread(&evt->energy, 2, 1, inFile); fread(&hit->energy, 2, 1, inFile);
fread(&evt->timestamp, 6, 1, inFile); if( hit->DPPType == DPPType::PSD ) fread(&hit->energy_short, 2, 1, inFile);
fread(&evt->fine_timestamp, 2, 1, inFile); fread(&hit->timestamp, 6, 1, inFile);
fread(&evt->flags_high_priority, 1, 1, inFile); fread(&hit->fine_timestamp, 2, 1, inFile);
fread(&evt->flags_low_priority, 2, 1, inFile); fread(&hit->flags_high_priority, 1, 1, inFile);
fread(&evt->downSampling, 1, 1, inFile); fread(&hit->flags_low_priority, 2, 1, inFile);
fread(&evt->board_fail, 1, 1, inFile); fread(&hit->downSampling, 1, 1, inFile);
fread(&evt->flush, 1, 1, inFile); fread(&hit->board_fail, 1, 1, inFile);
fread(&evt->trigger_threashold, 2, 1, inFile); fread(&hit->flush, 1, 1, inFile);
fread(&evt->event_size, 8, 1, inFile); fread(&hit->trigger_threashold, 2, 1, inFile);
fread(&evt->aggCounter, 4, 1, inFile); fread(&hit->event_size, 8, 1, inFile);
fread(&hit->aggCounter, 4, 1, inFile);
}else{ }else{
fseek(inFile, 31, SEEK_CUR); fseek(inFile, hit->DPPType == DPPType::PHA ? 31 : 33, SEEK_CUR);
} }
fread(&evt->traceLenght, 8, 1, inFile); fread(&hit->traceLenght, 8, 1, inFile);
if( isSkip == 0){ if( isSkip == 0){
fread(evt->analog_probes_type, 2, 1, inFile); fread(hit->analog_probes_type, 2, 1, inFile);
fread(evt->digital_probes_type, 4, 1, inFile); fread(hit->digital_probes_type, 4, 1, inFile);
fread(evt->analog_probes[0], evt->traceLenght*4, 1, inFile); fread(hit->analog_probes[0], hit->traceLenght*4, 1, inFile);
fread(evt->analog_probes[1], evt->traceLenght*4, 1, inFile); fread(hit->analog_probes[1], hit->traceLenght*4, 1, inFile);
fread(evt->digital_probes[0], evt->traceLenght, 1, inFile); fread(hit->digital_probes[0], hit->traceLenght, 1, inFile);
fread(evt->digital_probes[1], evt->traceLenght, 1, inFile); fread(hit->digital_probes[1], hit->traceLenght, 1, inFile);
fread(evt->digital_probes[2], evt->traceLenght, 1, inFile); fread(hit->digital_probes[2], hit->traceLenght, 1, inFile);
fread(evt->digital_probes[3], evt->traceLenght, 1, inFile); fread(hit->digital_probes[3], hit->traceLenght, 1, inFile);
}else{ }else{
fseek(inFile, 6 + evt->traceLenght*(12), SEEK_CUR); fseek(inFile, 6 + hit->traceLenght*(12), SEEK_CUR);
} }
}else if( evt->dataType == 1){ }else if( hit->dataType == DataFormat::OneTrace){
if( isSkip == 0 ){ if( isSkip == 0 ){
fread(&evt->channel, 1, 1, inFile); fread(&hit->channel, 1, 1, inFile);
fread(&evt->energy, 2, 1, inFile); fread(&hit->energy, 2, 1, inFile);
fread(&evt->timestamp, 6, 1, inFile); if( hit->DPPType == DPPType::PSD ) fread(&hit->energy_short, 2, 1, inFile);
fread(&evt->fine_timestamp, 2, 1, inFile); fread(&hit->timestamp, 6, 1, inFile);
fread(&evt->flags_high_priority, 1, 1, inFile); fread(&hit->fine_timestamp, 2, 1, inFile);
fread(&evt->flags_low_priority, 2, 1, inFile); fread(&hit->flags_high_priority, 1, 1, inFile);
fread(&hit->flags_low_priority, 2, 1, inFile);
}else{ }else{
fseek(inFile, 14, SEEK_CUR); fseek(inFile, hit->DPPType == DPPType::PHA ? 14 : 16, SEEK_CUR);
} }
fread(&evt->traceLenght, 8, 1, inFile); fread(&hit->traceLenght, 8, 1, inFile);
if( isSkip == 0){ if( isSkip == 0){
fread(&evt->analog_probes_type[0], 1, 1, inFile); fread(&hit->analog_probes_type[0], 1, 1, inFile);
fread(evt->analog_probes[0], evt->traceLenght*4, 1, inFile); fread(hit->analog_probes[0], hit->traceLenght*4, 1, inFile);
}else{ }else{
fseek(inFile, 1 + evt->traceLenght*4, SEEK_CUR); fseek(inFile, 1 + hit->traceLenght*4, SEEK_CUR);
} }
}else if( evt->dataType == 2){ }else if( hit->dataType == DataFormat::NoTrace){
if( isSkip == 0 ){ if( isSkip == 0 ){
fread(&evt->channel, 1, 1, inFile); fread(&hit->channel, 1, 1, inFile);
fread(&evt->energy, 2, 1, inFile); fread(&hit->energy, 2, 1, inFile);
fread(&evt->timestamp, 6, 1, inFile); if( hit->DPPType == DPPType::PSD ) fread(&hit->energy_short, 2, 1, inFile);
fread(&evt->fine_timestamp, 2, 1, inFile); fread(&hit->timestamp, 6, 1, inFile);
fread(&evt->flags_high_priority, 1, 1, inFile); fread(&hit->fine_timestamp, 2, 1, inFile);
fread(&evt->flags_low_priority, 2, 1, inFile); fread(&hit->flags_high_priority, 1, 1, inFile);
fread(&hit->flags_low_priority, 2, 1, inFile);
}else{ }else{
fseek(inFile, 14, SEEK_CUR); fseek(inFile, hit->DPPType == DPPType::PHA ? 14 : 16, SEEK_CUR);
} }
}else if( evt->dataType == 3){ }else if( hit->dataType == DataFormat::Minimum){
if( isSkip == 0 ){ if( isSkip == 0 ){
fread(&evt->channel, 1, 1, inFile); fread(&hit->channel, 1, 1, inFile);
fread(&evt->energy, 2, 1, inFile); fread(&hit->energy, 2, 1, inFile);
fread(&evt->timestamp, 6, 1, inFile); if( hit->DPPType == DPPType::PSD ) fread(&hit->energy_short, 2, 1, inFile);
fread(&hit->timestamp, 6, 1, inFile);
}else{ }else{
fseek(inFile, 9, SEEK_CUR); fseek(inFile, hit->DPPType == DPPType::PHA ? 9 : 11, SEEK_CUR);
} }
}else if( evt->dataType == 15){ }else if( hit->dataType == DataFormat::RAW){
fread(&evt->dataSize, 8, 1, inFile); fread(&hit->dataSize, 8, 1, inFile);
if( isSkip == 0){ if( isSkip == 0){
fread(evt->data, evt->dataSize, 1, inFile); fread(hit->data, hit->dataSize, 1, inFile);
}else{ }else{
fseek(inFile, evt->dataSize, SEEK_CUR); fseek(inFile, hit->dataSize, SEEK_CUR);
} }
} }
blockID ++; numBlock ++;
filePos = ftell(inFile); filePos = ftell(inFile);
return 0; return 0;
} }
void SolReader::RewindFile(){ void SolReader::RewindFile(){
rewind(inFile); rewind(inFile);
filePos = 0; filePos = 0;
blockID = 0; numBlock = 0;
} }
void SolReader::ScanNumBlock(){ void SolReader::ScanNumBlock(){
if( inFile == NULL ) return; if( inFile == NULL ) return;
if( feof(inFile) ) return; if( feof(inFile) ) return;
blockID = -1; numBlock = 0;
blockPos.clear(); blockPos.clear();
blockPos.push_back(0); blockPos.push_back(0);
totNumBlock = 0;
while( ReadNextBlock(1) == 0){ while( ReadNextBlock(1) == 0){
blockPos.push_back(filePos); blockPos.push_back(filePos);
printf("%ld, %.2f%% %u/%u\n\033[A\r", blockID, filePos*100./inFileSize, filePos, inFileSize); printf("%u, %.2f%% %u/%u\n\033[A\r", numBlock, filePos*100./inFileSize, filePos, inFileSize);
totNumBlock ++;
} }
blockID = -1; totNumBlock = numBlock;
numBlock = 0;
isScanned = true; isScanned = true;
printf("\nScan complete: number of data Block : %u\n", totNumBlock); printf("\nScan complete: number of data Block : %u\n", totNumBlock);
rewind(inFile); rewind(inFile);

169
working/script_multi.C Normal file
View File

@ -0,0 +1,169 @@
#include "../armory/SolReader.h"
#include "TH1.h"
#include "TMath.h"
#include "TH2.h"
#include "TStyle.h"
#include "TCanvas.h"
#include "TGraph.h"
void script_multi(std::string run = "033"){
SolReader* reader0 = new SolReader("../data_raw/*_" + run + "_00_21245_000.sol");
//SolReader* reader1 = new SolReader("/home/ryan/analysis/data_raw/test_" + run + "_01_21233_000.sol");
Hit * evt0 = reader0->hit;
//Event * evt1 = reader1->evt;
printf("----------file size: %u Byte\n", reader0->GetFileSize());
//printf("----------file size: %u Byte\n", reader1->GetFileSize());
reader0->ScanNumBlock();
//reader1->ScanNumBlock();
unsigned long startTime, endTime;
if( reader0->GetTotalNumBlock() == 0 ) return;
reader0->ReadBlock(0);
startTime = evt0->timestamp;
reader0->ReadBlock(reader0->GetTotalNumBlock() - 1);
endTime = evt0->timestamp;
double duration = double(endTime - startTime)*8./1e9;
printf("============== %lu ns = %.4f sec.\n", (endTime - startTime)*8, duration);
reader0->RewindFile();
// if( reader1->GetTotalNumBlock() == 0 ) return;
// reader1->ReadBlock(0);
// startTime = evt1->timestamp;
// reader1->ReadBlock(reader1->GetTotalNumBlock() - 1);
// endTime = evt1->timestamp;
// duration = double(endTime - startTime)*8./1e9;
// printf("============== %lu ns = %.4f sec.\n", (endTime - startTime)*8, duration);
// reader1->RewindFile();
//int minBlock = std::min(reader0->GetTotalNumBlock(), reader1->GetTotalNumBlock());
int minBlock = reader0->GetTotalNumBlock();
TH1F * hID = new TH1F("hID", "hID", 64, 0, 64);
TH1F * hTdiff = new TH1F("hTdiff", "tdiff", 300, -100, 200);
TH2F * hTdiff2 = new TH2F("hTdiff2", "tdiff vs evt", 400, 0, minBlock, 100, 0, 100);
TH1F * hRate0 = new TH1F("hRate0", "Rate", 20, 0, 20);
TH1F * hRate1 = new TH1F("hRate1", "Rate", 20, 0, 20); hRate1->SetLineColor(2);
TH1F * hE0 = new TH1F("hE0", "Energy", 400, 0, 30000);
TH1F * hE1 = new TH1F("hE1", "Energy", 400, 0, 30000); hE1->SetLineColor(2);
TH1F *hMulti = new TH1F("hMulti", "Multiplicy", 10, 0, 10);
std::vector<std::pair<int, uint64_t>> ts ;
for( int i = 0; i < minBlock; i++){
reader0->ReadNextBlock();
//reader1->ReadNextBlock();
if( i < 10 ){
printf("#################################################\n");
evt0->PrintAll();
}
ts.push_back(std::pair(evt0->channel,evt0->timestamp));
// printf("---------------------------------------------\n");
// evt1->PrintAll();
//if( evt0->channel == 30 ) evt0->PrintAll();
hID->Fill(evt0->channel);
if( evt0->channel == 0 ) {
hRate0->Fill( evt0->timestamp * 8 / 1e9);
hE0->Fill(evt0->energy);
}
if( evt0->channel == 30 ) {
hRate1->Fill( evt0->timestamp * 8 / 1e9);
hE1->Fill(evt0->energy);
}
//if( i < 10 ) printf("t0 : %10lu, t1 : %10lu, %10lu \n", evt0->timestamp, evt1->timestamp,
// evt0->timestamp > evt1->timestamp ? evt0->timestamp - evt1->timestamp : evt1->timestamp - evt0->timestamp);
//
//hTdiff->Fill(evt0->timestamp > evt1->timestamp ? evt0->timestamp - evt1->timestamp : evt1->timestamp - evt0->timestamp);
//hTdiff2->Fill(i, evt0->timestamp > evt1->timestamp ? evt0->timestamp - evt1->timestamp : evt1->timestamp - evt0->timestamp);
//if( i > 10 ) break;
}
delete reader0;
//delete reader1;
//build event
int coinWin = 200;
std::vector<std::vector<std::pair<int, uint64_t>>> events;
std::vector<int> multi;
events.push_back({ts[0]});
multi.push_back(1);
// Iterate through the vector starting from the second pair
for (size_t i = 1; i < ts.size(); ++i) {
uint64_t currentTimestamp = ts[i].second;
uint64_t previousTimestamp = ts[i - 1].second;
// Check if the timestamp difference is within coinWin
if (currentTimestamp - previousTimestamp <= coinWin) {
events.back().push_back(ts[i]); // Add to the current event
multi.back() += 1;
} else {
// Create a new event
events.push_back({ts[i]});
multi.push_back(1);
}
}
// Print the events
double maxTD = -999, minTD = 999;
for (size_t i = 0; i < events.size(); ++i) {
hMulti->Fill(multi[i]);
if( multi[i] > 1 ) {
printf("Event %zu, Multi : %d\n", i, multi[i]);
double haha[2];
for (size_t j = 0; j < events[i].size(); ++j) {
printf("Channel: %2d, Timestamp: %lu\n", events[i][j].first, events[i][j].second);
if( events[i][j].first == 0 ) haha[0] = static_cast<double> ( events[i][j].second );
if( events[i][j].first == 30 ) haha[1] = static_cast<double> ( events[i][j].second );
}
double TD = haha[1]-haha[0];
if(TD > maxTD && TD < 105) maxTD = TD;
if( TD < minTD && TD > -105) minTD = TD;
hTdiff->Fill( TD );
printf("-------------------\n");
}
}
printf(" max TD : %f \n", maxTD);
printf(" min TD : %f \n", minTD);
printf(" spn TD : %f \n", maxTD - minTD);
TCanvas * canvas = new TCanvas("canvas", "canvas", 1200, 800);
gStyle->SetOptStat("neiou");
canvas->Divide(3,2);
canvas->cd(1);hID->Draw();
canvas->cd(2);hE0->Draw(); hE1->Draw("same");
canvas->cd(3);hRate0->Draw(); hRate1->Draw("same");
canvas->cd(4);hTdiff->Draw();
canvas->cd(5); hMulti->Draw();
canvas->cd(5)->SetLogy(true);
canvas->cd(6); hE1->Draw();
//canvas->cd(2);hTdiff2->Draw();
}

49
working/script.C → working/script_single.C
View File

@ -3,10 +3,10 @@
#include "TH2.h" #include "TH2.h"
#include "TCanvas.h" #include "TCanvas.h"
void script(){ void script_single(std::string fileName){
SolReader * reader = new SolReader("../data_raw/Master_003_00_21245_000.sol"); SolReader * reader = new SolReader(fileName);
Event * evt = reader->evt; Hit * hit = reader->hit;
reader->ScanNumBlock(); reader->ScanNumBlock();
@ -14,27 +14,12 @@ void script(){
for( int i = 1; i < 10; i ++ ){ for( int i = 1; i < 10; i ++ ){
reader->ReadBlock(numBlock-i); reader->ReadBlock(numBlock-i);
evt->PrintEnergyTimeStamp(); hit->PrintEnergyTimeStamp();
} }
printf("======================== \n");
SolReader * reader1 = new SolReader("../data_raw/Master_003_01_21233_000.sol");
Event * evt1 = reader1->evt;
reader1->ScanNumBlock();
long numBlock1 = reader1->GetTotalNumBlock();
for( int i = 1; i < 10; i ++ ){
reader1->ReadBlock(numBlock1-i);
evt1->PrintEnergyTimeStamp();
}
/* /*
SolReader * reader = new SolReader("../data_raw/Master_003_00_21245_000.sol"); SolReader * reader = new SolReader("../data_raw/Master_003_00_21245_000.sol");
Event * evt = reader->evt; Event * hit = reader->hit;
printf("========== file size: %u Byte\n", reader->GetFileSize()); printf("========== file size: %u Byte\n", reader->GetFileSize());
@ -44,9 +29,9 @@ void script(){
unsigned long startTime, endTime; unsigned long startTime, endTime;
reader->ReadBlock(0); reader->ReadBlock(0);
startTime = evt->timestamp; startTime = hit->timestamp;
reader->ReadBlock(reader->GetTotalNumBlock() - 1); reader->ReadBlock(reader->GetTotalNumBlock() - 1);
endTime = evt->timestamp; endTime = hit->timestamp;
double duration = double(endTime - startTime)*8./1e9; double duration = double(endTime - startTime)*8./1e9;
printf("============== %lu ns = %.4f sec.\n", (endTime - startTime)*8, duration); printf("============== %lu ns = %.4f sec.\n", (endTime - startTime)*8, duration);
@ -70,16 +55,16 @@ void script(){
printf("########################## nBlock : %u, %u/%u\n", reader->GetBlockID(), printf("########################## nBlock : %u, %u/%u\n", reader->GetBlockID(),
reader->GetFilePos(), reader->GetFilePos(),
reader->GetFileSize()); reader->GetFileSize());
evt->PrintAll(); hit->PrintAll();
//evt->PrintAllTrace(); //hit->PrintAllTrace();
} }
h1->Fill(evt->timestamp); h1->Fill(hit->timestamp);
h2->Fill(evt->timestamp, i); h2->Fill(hit->timestamp, i);
if( i > 0 ){ if( i > 0 ){
if( evt->timestamp < tOld) printf("-------- time not sorted."); if( hit->timestamp < tOld) printf("-------- time not sorted.");
tOld = evt->timestamp; tOld = hit->timestamp;
} }
} }
@ -88,17 +73,17 @@ void script(){
h2->Draw(); h2->Draw();
*/ */
//printf("reader traceLength : %lu \n", evt->traceLenght); //printf("reader traceLength : %lu \n", hit->traceLenght);
/* /*
for( int i = 0; i < evt->traceLenght; i++){ for( int i = 0; i < hit->traceLenght; i++){
printf("%4d| %d\n", i, evt->analog_probes[0][i]); printf("%4d| %d\n", i, hit->analog_probes[0][i]);
} }
*/ */
evt = NULL; hit = NULL;
delete reader; delete reader;
} }