diff --git a/.gitignore b/.gitignore index 98b4f2b..6856f48 100644 --- a/.gitignore +++ b/.gitignore @@ -6,6 +6,8 @@ EventBuilder *.d *.so +.gdb_history + data_raw root_data diff --git a/.vscode/settings.json b/.vscode/settings.json index 1e65214..6b63c3c 100644 --- a/.vscode/settings.json +++ b/.vscode/settings.json @@ -119,7 +119,9 @@ "locale": "cpp", "__verbose_abort": "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, diff --git a/README.md b/README.md new file mode 100644 index 0000000..8950318 --- /dev/null +++ b/README.md @@ -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. + + diff --git a/armory/EventBuilder.cpp b/armory/EventBuilder.cpp index 391d4e8..18244be 100644 --- a/armory/EventBuilder.cpp +++ b/armory/EventBuilder.cpp @@ -10,11 +10,11 @@ //#include "TClonesArray.h" // plan to save trace as TVector with TClonesArray //#include "TVector.h" -#define MAX_ID 64 +#define MAX_MULTI 64 #define MAX_TRACE_LEN 2500 SolReader ** reader; -Event ** evt; +Hit ** hit; unsigned long totFileSize = 0; unsigned long processedFileSize = 0; @@ -25,48 +25,51 @@ std::vector> group; // group[i][j], i = group ID, j = group mem void findEarliestTime(int &fileID, int & digiID){ - unsigned long firstTime = 0; + unsigned long firstTime = 0; for( int i = 0; i < (int) activeFileID.size(); i++){ int id = activeFileID[i]; if( i == 0 ) { - firstTime = evt[id]->timestamp; + firstTime = hit[id]->timestamp; fileID = id; - digiID = i; - //printf("%d | %ld %lu %d | %d \n", id, reader[id]->GetBlockID(), evt[id]->timestamp, evt[id]->channel, (int) activeFileID.size()); + //printf("%d | %d %lu %d | %d \n", id, reader[id]->GetBlockID(), hit[id]->timestamp, hit[id]->channel, (int) activeFileID.size()); continue; } - if( evt[id]->timestamp <= firstTime) { - firstTime = evt[id]->timestamp; + if( hit[id]->timestamp <= firstTime) { + firstTime = hit[id]->timestamp; fileID = id; - digiID = i; - //printf("%d | %ld %lu %d | %d \n", id, reader[id]->GetBlockID(), evt[id]->timestamp, evt[id]->channel, (int) activeFileID.size()); + //printf("%d | %d %lu %d | %d \n", id, reader[id]->GetBlockID(), hit[id]->timestamp, hit[id]->channel, (int) activeFileID.size()); } } + } -unsigned long long evID = 0; -int multi = 0; -int bd[MAX_ID] = {0}; -int ch[MAX_ID] = {0}; -int e[MAX_ID] = {0}; -unsigned long long e_t[MAX_ID] = {0}; -unsigned short lowFlag[MAX_ID] = {0}; -unsigned short highFlag[MAX_ID] = {0}; -int traceLen[MAX_ID] = {0}; -int trace[MAX_ID][MAX_TRACE_LEN] = {0}; +unsigned long long evID = 0; +unsigned int multi = 0; +unsigned short bd[MAX_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}; //for PSD energy short +unsigned long long e_t[MAX_MULTI] = {0}; +unsigned short lowFlag[MAX_MULTI] = {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){ - bd[multi] = digiID; - ch[multi] = evt[fileID]->channel; - e[multi] = evt[fileID]->energy; - e_t[multi] = evt[fileID]->timestamp; - lowFlag[multi] = evt[fileID]->flags_low_priority; - highFlag[multi] = evt[fileID]->flags_high_priority; + bd[multi] = digiID; + sn[multi] = digiID; + ch[multi] = hit[fileID]->channel; + e[multi] = hit[fileID]->energy; + e2[multi] = hit[fileID]->energy_short; + e_t[multi] = hit[fileID]->timestamp; + lowFlag[multi] = hit[fileID]->flags_low_priority; + highFlag[multi] = hit[fileID]->flags_high_priority; 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++){ - 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 reader = new SolReader*[nFile]; - evt = new Event *[nFile]; + hit = new Hit *[nFile]; for( int i = 0 ; i < nFile ; i++){ reader[i] = new SolReader(inFileName[i].Data()); - evt[i] = reader[i]->evt; //TODO check is file open propertly - - //reader[i]->ScanNumBlock(); - + hit[i] = reader[i]->hit; //TODO check is file open propertly reader[i]->ReadNextBlock(); // read the first block } @@ -144,16 +144,20 @@ int main(int argc, char ** argv){ int digiID = f1.Remove(pos).Atoi(); 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); pos = fn.First('.'); // get the file id; - int fileID = fn.Remove(pos).Atoi(); + int indexID = fn.Remove(pos).Atoi(); - std::vector haha = {i, digiID, fileID}; + int fileID = i; + std::vector haha = {fileID, digiID, indexID, digisn}; idList.push_back(haha); } + // sort by digiID std::sort(idList.begin(), idList.end(), [](const std::vector& a, const std::vector& b){ if (a[1] == b[1]) { return a[2] < b[2]; @@ -161,11 +165,11 @@ int main(int argc, char ** argv){ return a[1] < b[1]; }); - group.clear(); + group.clear(); // group[i][j], i is the group Index = digiID int last_id = 0; std::vector kaka; for( int i = 0; i < (int) idList.size() ; i++){ - if( i == 0 ) { + if( i == 0 ) { kaka.clear(); last_id = idList[i][1]; kaka.push_back(idList[i][0]); @@ -208,10 +212,12 @@ int main(int argc, char ** argv){ TTree * tree = new TTree("tree", outFileName); tree->Branch("evID", &evID, "event_ID/l"); - tree->Branch("multi", &multi, "multi/I"); - tree->Branch("bd", bd, "board[multi]/I"); - tree->Branch("ch", ch, "channel[multi]/I"); - tree->Branch("e", e, "energy[multi]/I"); + tree->Branch("multi", &multi, "multi/i"); + tree->Branch("bd", bd, "board[multi]/s"); + tree->Branch("sn", sn, "sn[multi]/s"); + 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("lowFlag", lowFlag, "lowFlag[multi]/s"); tree->Branch("highFlag", highFlag, "highFlag[multi]/s"); @@ -241,15 +247,16 @@ int main(int argc, char ** argv){ findEarliestTime(fileID, digiID); fillData(fileID, digiID, saveTrace); - unsigned long firstTimeStamp = evt[fileID]->timestamp; + unsigned long firstTimeStamp = hit[fileID]->timestamp; unsigned long lastTimeStamp = 0; - int last_precentage = 0; while((activeFileID.size() > 0)){ 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); }else{ outRootFile->cd(); @@ -259,7 +266,7 @@ int main(int argc, char ** argv){ multi = 0; fillData(fileID, digiID, saveTrace); } - + ///========= check is file finished. if( relay[digiID]){ groupIndex[digiID] ++; @@ -271,8 +278,7 @@ int main(int argc, char ** argv){ } relay[digiID] = false; } - if( reader[fileID]->IsEndOfFile() ) relay[digiID] = true; - + ///========= calculate progress processedFileSize = 0; for( int p = 0; p < (int) group.size(); p ++){ @@ -302,7 +308,7 @@ int main(int argc, char ** argv){ double percentage = processedFileSize * 100/ totFileSize; printf("Processed : %llu, %.0f%% | %lu/%lu \n", evID, percentage, processedFileSize, totFileSize); - lastTimeStamp = evt[fileID]->timestamp; + lastTimeStamp = hit[fileID]->timestamp; //*=========================================== save file outRootFile->cd(); tree->Fill(); diff --git a/armory/Event.h b/armory/Hit.h similarity index 59% rename from armory/Event.h rename to armory/Hit.h index 2a741c7..948a1de 100644 --- a/armory/Event.h +++ b/armory/Hit.h @@ -1,5 +1,5 @@ -#ifndef EVENT_H -#define EVENT_H +#ifndef HIT_H +#define HIT_H #include #include @@ -8,27 +8,46 @@ #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: unsigned short dataType; + std::string DPPType; ///============= for dpp-pha - uint8_t channel; // 6 bit - uint16_t energy; // 16 bit - uint64_t timestamp; // 48 bit + uint8_t channel; // 6 bit + uint16_t energy; // 16 bit + uint16_t energy_short; // 16 bit, only for PSD + uint64_t timestamp; // 48 bit uint16_t fine_timestamp; // 16 bit uint16_t flags_low_priority; // 12 bit uint16_t flags_high_priority; // 8 bit - size_t traceLenght; // 64 bit - uint8_t downSampling; // 8 bit + size_t traceLenght; // 64 bit + uint8_t downSampling; // 8 bit bool board_fail; bool flush; - uint8_t analog_probes_type[2]; // 3 bit - uint8_t digital_probes_type[4]; // 4 bit - int32_t * analog_probes[2]; // 18 bit - uint8_t * digital_probes[4]; // 1 bit - uint16_t trigger_threashold; // 16 bit + uint8_t analog_probes_type[2]; // 3 bit for PHA, 4 bit for PSD + uint8_t digital_probes_type[4]; // 4 bit for PHA, 5 bit for PSD + int32_t * analog_probes[2]; // 18 bit + uint8_t * digital_probes[4]; // 1 bit + uint16_t trigger_threashold; // 16 bit size_t event_size; // 64 bit uint32_t aggCounter; // 32 bit @@ -39,17 +58,21 @@ class Event { bool isTraceAllZero; - Event(){ + Hit(){ Init(); } - ~Event(){ + ~Hit(){ ClearMemory(); } void Init(){ + DPPType = DPPType::PHA; + dataType = DataFormat::ALL; + channel = 0; energy = 0; + energy_short = 0; timestamp = 0; fine_timestamp = 0; downSampling = 0; @@ -92,11 +115,12 @@ class Event { isTraceAllZero = true; } - void SetDataType(unsigned int type){ + void SetDataType(unsigned int type, std::string dppType){ dataType = type; + DPPType = dppType; ClearMemory(); - if( dataType == 0xF){ + if( dataType == DataFormat::RAW){ data = new uint8_t[20*1024*1024]; }else{ analog_probes[0] = new int32_t[MaxTraceLenght]; @@ -128,34 +152,68 @@ class Event { } 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){ - switch(probeType){ - case 0: return "ADC"; - case 1: return "Time filter"; - case 2: return "Energy filter"; - default : return "none"; + + if( DPPType == DPPType::PHA){ + switch(probeType){ + case 0: return "ADC"; + 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){ - switch(probeType){ - case 0: return "Trigger"; - case 1: return "Time filter armed"; - case 2: return "Re-trigger guard"; - case 3: return "Energy filter baseline freeze"; - case 4: return "Energy filter peaking"; - case 5: return "Energy filter peaking ready"; - case 6: return "Energy filter pile-up guard"; - case 7: return "Event pile-up"; - case 8: return "ADC saturation"; - case 9: return "ADC saturation protection"; - case 10: return "Post-saturation event"; - case 11: return "Energy filter saturation"; - case 12: return "Signal inhibit"; - default : return "none"; + + if( DPPType == DPPType::PHA){ + switch(probeType){ + case 0: return "Trigger"; + case 1: return "Time filter armed"; + case 2: return "Re-trigger guard"; + case 3: return "Energy filter baseline freeze"; + case 4: return "Energy filter peaking"; + case 5: return "Energy filter peaking ready"; + case 6: return "Energy filter pile-up guard"; + case 7: return "Event pile-up"; + case 8: return "ADC saturation"; + case 9: return "ADC saturation protection"; + case 10: return "Post-saturation event"; + 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 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("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("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], diff --git a/armory/Makefile b/armory/Makefile index ab42e06..4de3a66 100644 --- a/armory/Makefile +++ b/armory/Makefile @@ -1,10 +1,11 @@ CC=g++ +CFLAG= -g ROOTFLAG=`root-config --cflags --glibs` all: EventBuilder -EventBuilder: EventBuilder.cpp SolReader.h - $(CC) EventBuilder.cpp -o EventBuilder ${ROOTFLAG} +EventBuilder: EventBuilder.cpp SolReader.h Hit.h + $(CC) $(CFLAG) EventBuilder.cpp -o EventBuilder ${ROOTFLAG} clean: diff --git a/armory/SolReader.h b/armory/SolReader.h index 32f3c24..424405f 100644 --- a/armory/SolReader.h +++ b/armory/SolReader.h @@ -8,7 +8,7 @@ #include #include // 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 @@ -20,7 +20,7 @@ class SolReader { unsigned int totNumBlock; unsigned short blockStartIdentifier; - long blockID; + unsigned int numBlock; bool isScanned; void init(); @@ -29,7 +29,7 @@ class SolReader { public: 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(); void OpenFile(std::string fileName); @@ -38,25 +38,25 @@ class SolReader { 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 GetFilePos() const {return filePos;} unsigned int GetFileSize() const {return inFileSize;} - - bool IsEndOfFile() {return (filePos >= inFileSize ? true : false);} - + void RewindFile(); - Event * evt; + Hit * hit; }; void SolReader::init(){ inFileSize = 0; - blockID = -1; + numBlock = 0; filePos = 0; totNumBlock = 0; - evt = new Event(); + hit = new Hit(); isScanned = false; @@ -68,16 +68,15 @@ SolReader::SolReader(){ init(); } -SolReader::SolReader(std::string fileName, unsigned short dataType){ +SolReader::SolReader(std::string fileName, unsigned short dataType = 0){ init(); OpenFile(fileName); - evt->SetDataType(dataType); + hit->SetDataType(dataType, DPPType::PHA); } SolReader::~SolReader(){ - //printf("%s\n", __func__); if( !inFile ) fclose(inFile); - delete evt; + delete hit; } 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]); fseek(inFile, blockPos[index], SEEK_CUR); + filePos = blockPos[index]; - blockID = index; - blockID --; + numBlock = index; + return ReadNextBlock(); } inline int SolReader::ReadNextBlock(int isSkip){ if( inFile == NULL ) return -1; - if( feof(inFile) ) return -1; + if( feof(inFile) ) return -1; if( filePos >= inFileSize) return -1; fread(&blockStartIdentifier, 2, 1, inFile); - if( (blockStartIdentifier & 0xAAA0) != 0xAAA0 ) { + if( (blockStartIdentifier & 0xAA00) != 0xAA00 ) { printf("header fail.\n"); return -2 ; } - if( ( blockStartIdentifier & 0xF ) == 15 ){ - evt->SetDataType(15); + if( ( blockStartIdentifier & 0xF ) == DataFormat::RAW ){ + 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 ){ - fread(&evt->channel, 1, 1, inFile); - fread(&evt->energy, 2, 1, inFile); - fread(&evt->timestamp, 6, 1, inFile); - fread(&evt->fine_timestamp, 2, 1, inFile); - fread(&evt->flags_high_priority, 1, 1, inFile); - fread(&evt->flags_low_priority, 2, 1, inFile); - fread(&evt->downSampling, 1, 1, inFile); - fread(&evt->board_fail, 1, 1, inFile); - fread(&evt->flush, 1, 1, inFile); - fread(&evt->trigger_threashold, 2, 1, inFile); - fread(&evt->event_size, 8, 1, inFile); - fread(&evt->aggCounter, 4, 1, inFile); + fread(&hit->channel, 1, 1, inFile); + fread(&hit->energy, 2, 1, inFile); + if( hit->DPPType == DPPType::PSD ) fread(&hit->energy_short, 2, 1, inFile); + fread(&hit->timestamp, 6, 1, inFile); + fread(&hit->fine_timestamp, 2, 1, inFile); + fread(&hit->flags_high_priority, 1, 1, inFile); + fread(&hit->flags_low_priority, 2, 1, inFile); + fread(&hit->downSampling, 1, 1, inFile); + fread(&hit->board_fail, 1, 1, inFile); + fread(&hit->flush, 1, 1, inFile); + fread(&hit->trigger_threashold, 2, 1, inFile); + fread(&hit->event_size, 8, 1, inFile); + fread(&hit->aggCounter, 4, 1, inFile); }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){ - fread(evt->analog_probes_type, 2, 1, inFile); - fread(evt->digital_probes_type, 4, 1, inFile); - fread(evt->analog_probes[0], evt->traceLenght*4, 1, inFile); - fread(evt->analog_probes[1], evt->traceLenght*4, 1, inFile); - fread(evt->digital_probes[0], evt->traceLenght, 1, inFile); - fread(evt->digital_probes[1], evt->traceLenght, 1, inFile); - fread(evt->digital_probes[2], evt->traceLenght, 1, inFile); - fread(evt->digital_probes[3], evt->traceLenght, 1, inFile); + fread(hit->analog_probes_type, 2, 1, inFile); + fread(hit->digital_probes_type, 4, 1, inFile); + fread(hit->analog_probes[0], hit->traceLenght*4, 1, inFile); + fread(hit->analog_probes[1], hit->traceLenght*4, 1, inFile); + fread(hit->digital_probes[0], hit->traceLenght, 1, inFile); + fread(hit->digital_probes[1], hit->traceLenght, 1, inFile); + fread(hit->digital_probes[2], hit->traceLenght, 1, inFile); + fread(hit->digital_probes[3], hit->traceLenght, 1, inFile); }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 ){ - fread(&evt->channel, 1, 1, inFile); - fread(&evt->energy, 2, 1, inFile); - fread(&evt->timestamp, 6, 1, inFile); - fread(&evt->fine_timestamp, 2, 1, inFile); - fread(&evt->flags_high_priority, 1, 1, inFile); - fread(&evt->flags_low_priority, 2, 1, inFile); + fread(&hit->channel, 1, 1, inFile); + fread(&hit->energy, 2, 1, inFile); + if( hit->DPPType == DPPType::PSD ) fread(&hit->energy_short, 2, 1, inFile); + fread(&hit->timestamp, 6, 1, inFile); + fread(&hit->fine_timestamp, 2, 1, inFile); + fread(&hit->flags_high_priority, 1, 1, inFile); + fread(&hit->flags_low_priority, 2, 1, inFile); }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){ - fread(&evt->analog_probes_type[0], 1, 1, inFile); - fread(evt->analog_probes[0], evt->traceLenght*4, 1, inFile); + fread(&hit->analog_probes_type[0], 1, 1, inFile); + fread(hit->analog_probes[0], hit->traceLenght*4, 1, inFile); }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 ){ - fread(&evt->channel, 1, 1, inFile); - fread(&evt->energy, 2, 1, inFile); - fread(&evt->timestamp, 6, 1, inFile); - fread(&evt->fine_timestamp, 2, 1, inFile); - fread(&evt->flags_high_priority, 1, 1, inFile); - fread(&evt->flags_low_priority, 2, 1, inFile); + fread(&hit->channel, 1, 1, inFile); + fread(&hit->energy, 2, 1, inFile); + if( hit->DPPType == DPPType::PSD ) fread(&hit->energy_short, 2, 1, inFile); + fread(&hit->timestamp, 6, 1, inFile); + fread(&hit->fine_timestamp, 2, 1, inFile); + fread(&hit->flags_high_priority, 1, 1, inFile); + fread(&hit->flags_low_priority, 2, 1, inFile); }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 ){ - fread(&evt->channel, 1, 1, inFile); - fread(&evt->energy, 2, 1, inFile); - fread(&evt->timestamp, 6, 1, inFile); + fread(&hit->channel, 1, 1, inFile); + fread(&hit->energy, 2, 1, inFile); + if( hit->DPPType == DPPType::PSD ) fread(&hit->energy_short, 2, 1, inFile); + fread(&hit->timestamp, 6, 1, inFile); }else{ - fseek(inFile, 9, SEEK_CUR); + fseek(inFile, hit->DPPType == DPPType::PHA ? 9 : 11, SEEK_CUR); } - }else if( evt->dataType == 15){ - fread(&evt->dataSize, 8, 1, inFile); + }else if( hit->dataType == DataFormat::RAW){ + fread(&hit->dataSize, 8, 1, inFile); if( isSkip == 0){ - fread(evt->data, evt->dataSize, 1, inFile); + fread(hit->data, hit->dataSize, 1, inFile); }else{ - fseek(inFile, evt->dataSize, SEEK_CUR); + fseek(inFile, hit->dataSize, SEEK_CUR); } } - blockID ++; + numBlock ++; filePos = ftell(inFile); - return 0; } void SolReader::RewindFile(){ rewind(inFile); filePos = 0; - blockID = 0; + numBlock = 0; } void SolReader::ScanNumBlock(){ if( inFile == NULL ) return; if( feof(inFile) ) return; - blockID = -1; + numBlock = 0; blockPos.clear(); blockPos.push_back(0); - totNumBlock = 0; while( ReadNextBlock(1) == 0){ blockPos.push_back(filePos); - printf("%ld, %.2f%% %u/%u\n\033[A\r", blockID, filePos*100./inFileSize, filePos, inFileSize); - totNumBlock ++; + printf("%u, %.2f%% %u/%u\n\033[A\r", numBlock, filePos*100./inFileSize, filePos, inFileSize); } - blockID = -1; + totNumBlock = numBlock; + numBlock = 0; isScanned = true; printf("\nScan complete: number of data Block : %u\n", totNumBlock); rewind(inFile); diff --git a/working/script_multi.C b/working/script_multi.C new file mode 100644 index 0000000..37aa464 --- /dev/null +++ b/working/script_multi.C @@ -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> 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>> events; + std::vector 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 ( events[i][j].second ); + if( events[i][j].first == 30 ) haha[1] = static_cast ( 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(); + + +} \ No newline at end of file diff --git a/working/script.C b/working/script_single.C similarity index 62% rename from working/script.C rename to working/script_single.C index cc0f65e..99c9c3a 100644 --- a/working/script.C +++ b/working/script_single.C @@ -3,10 +3,10 @@ #include "TH2.h" #include "TCanvas.h" -void script(){ +void script_single(std::string fileName){ - SolReader * reader = new SolReader("../data_raw/Master_003_00_21245_000.sol"); - Event * evt = reader->evt; + SolReader * reader = new SolReader(fileName); + Hit * hit = reader->hit; reader->ScanNumBlock(); @@ -14,27 +14,12 @@ void script(){ for( int i = 1; i < 10; 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"); - Event * evt = reader->evt; + Event * hit = reader->hit; printf("========== file size: %u Byte\n", reader->GetFileSize()); @@ -44,9 +29,9 @@ void script(){ unsigned long startTime, endTime; reader->ReadBlock(0); - startTime = evt->timestamp; + startTime = hit->timestamp; reader->ReadBlock(reader->GetTotalNumBlock() - 1); - endTime = evt->timestamp; + endTime = hit->timestamp; double duration = double(endTime - startTime)*8./1e9; printf("============== %lu ns = %.4f sec.\n", (endTime - startTime)*8, duration); @@ -70,16 +55,16 @@ void script(){ printf("########################## nBlock : %u, %u/%u\n", reader->GetBlockID(), reader->GetFilePos(), reader->GetFileSize()); - evt->PrintAll(); - //evt->PrintAllTrace(); + hit->PrintAll(); + //hit->PrintAllTrace(); } - h1->Fill(evt->timestamp); - h2->Fill(evt->timestamp, i); + h1->Fill(hit->timestamp); + h2->Fill(hit->timestamp, i); if( i > 0 ){ - if( evt->timestamp < tOld) printf("-------- time not sorted."); - tOld = evt->timestamp; + if( hit->timestamp < tOld) printf("-------- time not sorted."); + tOld = hit->timestamp; } } @@ -88,17 +73,17 @@ void script(){ 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; } \ No newline at end of file