MLXID/DataProcess
0
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Add hists and logics for 3-photon clusters

Browse Source
This commit is contained in:
xie_x1
2026-05-29 16:06:11 +02:00
parent d4f5c7c5da
commit d9af7044ef
1 changed files with 95 additions and 1 deletions
Download Patch File Download Diff File Expand all files Collapse all files
UsefulFuncs.py
+95 -1
View File
@@ -12,6 +12,7 @@ _aduToKev3DMap = None ### global caliFile
_pedestalAduFrame = None ### global pedestalAduFrame
_noiseEneFrame = None ### global noiseEneFrame
nChunks = 16
clusterSize3Photon = 9
clusterSize2Photon = 7
clusterSize1Photon = 3
@@ -111,26 +112,51 @@ def bookHistograms(energy, suffix='', energyBinWidth=0.1, isMC=False):
)
histograms['h1_2PhotonClusterEnergy'].SetTitle('2-Photon Cluster Energy;Energy [keV];Counts')
histograms['h1_3PhotonClusterEnergy'] = TH1D(
f'h1_3PhotonClusterEnergy{suffix}', f'h1_3PhotonClusterEnergy{suffix}', nEnergyBins, -1, energy*6
)
histograms['h1_3PhotonClusterEnergy'].SetTitle('3-Photon Cluster Energy;Energy [keV];Counts')
histograms['h1_2PhotonClusterSize'] = TH1D(
f'h1_2PhotonClusterSize{suffix}', f'h1_2PhotonClusterSize{suffix}', 30, 0, 30
)
histograms['h1_2PhotonClusterSize'].SetTitle('2-Photon Cluster Size;Cluster Size [pixel];Counts')
histograms['h1_3PhotonClusterSize'] = TH1D(
f'h1_3PhotonClusterSize{suffix}', f'h1_3PhotonClusterSize{suffix}', 30, 0, 30
)
histograms['h1_3PhotonClusterSize'].SetTitle('3-Photon Cluster Size;Cluster Size [pixel];Counts')
histograms['h1_2PhotonClusterDiameterX'] = TH1D(
f'h1_2PhotonClusterDiameterX{suffix}', f'h1_2PhotonClusterDiameterX{suffix}', 10, 0, 10
)
histograms['h1_2PhotonClusterDiameterX'].SetTitle('2-Photon Cluster Diameter X;Cluster Diameter X [pixel];Counts')
histograms['h1_3PhotonClusterDiameterX'] = TH1D(
f'h1_3PhotonClusterDiameterX{suffix}', f'h1_3PhotonClusterDiameterX{suffix}', 15, 0, 15
)
histograms['h1_3PhotonClusterDiameterX'].SetTitle('3-Photon Cluster Diameter X;Cluster Diameter X [pixel];Counts')
histograms['h1_2PhotonClusterDiameterY'] = TH1D(
f'h1_2PhotonClusterDiameterY{suffix}', f'h1_2PhotonClusterDiameterY{suffix}', 10, 0, 10
)
histograms['h1_2PhotonClusterDiameterY'].SetTitle('2-Photon Cluster Diameter Y;Cluster Diameter Y [pixel];Counts')
histograms['h1_3PhotonClusterDiameterY'] = TH1D(
f'h1_3PhotonClusterDiameterY{suffix}', f'h1_3PhotonClusterDiameterY{suffix}', 15, 0, 15
)
histograms['h1_3PhotonClusterDiameterY'].SetTitle('3-Photon Cluster Diameter Y;Cluster Diameter Y [pixel];Counts')
histograms['h1_2PhotonClusterPixelEnergy'] = TH1D(
f'h1_2PhotonClusterPixelEnergy{suffix}', f'h1_2PhotonClusterPixelEnergy{suffix}', nEnergyBins, -1, energy*2.5
)
histograms['h1_2PhotonClusterPixelEnergy'].SetTitle('2-Photon Cluster Pixel Energy;Energy [keV];Counts')
histograms['h1_3PhotonClusterPixelEnergy'] = TH1D(
f'h1_3PhotonClusterPixelEnergy{suffix}', f'h1_3PhotonClusterPixelEnergy{suffix}', nEnergyBins, -1, energy*2.5
)
histograms['h1_3PhotonClusterPixelEnergy'].SetTitle('3-Photon Cluster Pixel Energy;Energy [keV];Counts')
histograms['h2_SumFrame'] = TH2D(
f'h2_SumFrame{suffix}', f'h2_SumFrame{suffix}', roi_width, roi_x0, roi_x1, roi_height, roi_y0, roi_y1
)
@@ -156,18 +182,28 @@ def bookHistograms(energy, suffix='', energyBinWidth=0.1, isMC=False):
)
histograms['h2_2PhotonSumSample'].SetTitle('2-Photon Sum Sample;Sum Sample;Counts')
histograms['h2_3PhotonSumSample'] = TH2D(
f'h2_3PhotonSumSample{suffix}', f'h2_3PhotonSumSample{suffix}', clusterSize3Photon, 0, clusterSize3Photon, clusterSize3Photon, 0, clusterSize3Photon
)
histograms['h2_3PhotonSumSample'].SetTitle('3-Photon Sum Sample;Sum Sample;Counts')
histograms['h2_2PhotonSumFrame'] = TH2D(
f'h2_2PhotonSumFrame{suffix}', f'h2_2PhotonSumFrame{suffix}', roi_width, roi_x0, roi_x1, roi_height, roi_y0, roi_y1
)
histograms['h2_2PhotonSumFrame'].SetTitle('2-Photon Hits Sum Frame;X [pixel];Y [pixel];Energy Sum [keV]')
histograms['h2_3PhotonSumFrame'] = TH2D(
f'h2_3PhotonSumFrame{suffix}', f'h2_3PhotonSumFrame{suffix}', roi_width, roi_x0, roi_x1, roi_height, roi_y0, roi_y1
)
histograms['h2_3PhotonSumFrame'].SetTitle('3-Photon Hits Sum Frame;X [pixel];Y [pixel];Energy Sum [keV]')
histograms['h2_ClusterSizeVsEnergy'] = TH2D(
f'h2_ClusterSizeVsEnergy{suffix}', f'h2_ClusterSizeVsEnergy{suffix}', nEnergyBins, -1, histMaxEnergy, 30, 0, 30
)
histograms['h2_ClusterSizeVsEnergy'].SetTitle('Cluster Size vs Energy;Energy [keV];Cluster Size [pixel]')
histograms['h2_DiameterVsEnergy'] = TH2D(
f'h2_DiameterVsEnergy{suffix}', f'h2_DiameterVsEnergy{suffix}', nEnergyBins, -1, histMaxEnergy, 10, 0, 10
f'h2_DiameterVsEnergy{suffix}', f'h2_DiameterVsEnergy{suffix}', nEnergyBins, -1, histMaxEnergy, 15, 0, 15
)
histograms['h2_DiameterVsEnergy'].SetTitle('Cluster Diameter vs Energy;Energy [keV];Cluster Diameter [pixel]')
@@ -224,6 +260,19 @@ def _processFrames(idxChunk): ### for both single and double photon events, usin
)
cluster_2Photon_list = []
refpoint_2Photon_list = []
h5_3Photon_file = h5py.File(f'{_cfg["outputFolder"]}/3Photon_CS{clusterSize3Photon}_chunk{idxChunk}.h5', 'w')
dset_3Photon_clusters = h5_3Photon_file.create_dataset(
'clusters', (0, clusterSize3Photon, clusterSize3Photon), maxshape=(None, clusterSize3Photon, clusterSize3Photon), dtype='f4',
chunks=True, compression='gzip'
)
dset_3Photon_refs = h5_3Photon_file.create_dataset(
'referencePoint', (0, 2), maxshape=(None, 2), dtype='i4',
chunks=True
)
cluster_3Photon_list = []
refpoint_3Photon_list = []
BUFFER_THRESHOLD = 10000
def flush_h5_buffer():
if not cluster_1Photon_list: return # no new clusters to write
@@ -236,6 +285,7 @@ def _processFrames(idxChunk): ### for both single and double photon events, usin
cluster_1Photon_list.clear()
refpoint_1Photon_list.clear()
if not cluster_2Photon_list: return # no new clusters to write
old_len = dset_2Photon_clusters.shape[0]
new_len = old_len + len(cluster_2Photon_list)
dset_2Photon_clusters.resize((new_len, clusterSize2Photon, clusterSize2Photon))
@@ -245,6 +295,16 @@ def _processFrames(idxChunk): ### for both single and double photon events, usin
cluster_2Photon_list.clear()
refpoint_2Photon_list.clear()
if not cluster_3Photon_list: return # no new clusters to write
old_len = dset_3Photon_clusters.shape[0]
new_len = old_len + len(cluster_3Photon_list)
dset_3Photon_clusters.resize((new_len, clusterSize3Photon, clusterSize3Photon))
dset_3Photon_refs.resize((new_len, 2))
dset_3Photon_clusters[old_len:new_len] = cluster_3Photon_list
dset_3Photon_refs[old_len:new_len] = refpoint_3Photon_list
cluster_3Photon_list.clear()
refpoint_3Photon_list.clear()
for idxFrame in range(startFrame, endFrame):
idxFile, idxFrame = divmod(idxFrame, nFramePerFile)
try:
@@ -357,6 +417,40 @@ def _processFrames(idxChunk): ### for both single and double photon events, usin
cluster_2Photon_list.append(cluster_2Photon)
refpoint_2Photon_list.append([int(ref_x)+Roi[0], int(ref_y)+Roi[2]])
if 3 * Energy - 3 * selectionRange < energy < 3 * Energy + 3 * selectionRange:
x_center = np.sum(xs * enes) / np.sum(enes)
y_center = np.sum(ys * enes) / np.sum(enes)
ref_x = int(x_center - clusterSize3Photon / 2) + 1 ### refered to the lower-left corner of the cluster
ref_y = int(y_center - clusterSize3Photon / 2) + 1
if int(ref_x) < 0 or int(ref_x)+clusterSize3Photon > signalEneFrame.shape[1] or int(ref_y) < 0 or int(ref_y)+clusterSize3Photon > signalEneFrame.shape[0]:
continue ### skip clusters too close to the border
for i in range(len(xs)):
_hists['h2_3PhotonSumFrame'].Fill(xs[i]+Roi[0], ys[i]+Roi[2], enes[i])
_hists['h1_3PhotonClusterPixelEnergy'].Fill(enes[i])
_hists['h1_3PhotonClusterSize'].Fill(clusterSize)
diameterX = max(xs) - min(xs) + 1
diameterY = max(ys) - min(ys) + 1
_hists['h1_3PhotonClusterDiameterX'].Fill(diameterX)
_hists['h1_3PhotonClusterDiameterY'].Fill(diameterY)
### filling the 3-photon cluster. Add background as well.
cluster_3Photon = np.zeros((clusterSize3Photon, clusterSize3Photon), dtype=np.float32)
for i in range(len(xs)):
x_rel = xs[i] - int(ref_x)
y_rel = ys[i] - int(ref_y)
if 0 <= x_rel < clusterSize3Photon and 0 <= y_rel < clusterSize3Photon:
cluster_3Photon[y_rel, x_rel] = enes[i]
_hists['h2_3PhotonSumSample'].Fill(x_rel, y_rel, enes[i])
cluster_3Photon += signalEneFrame[int(ref_y):int(ref_y)+clusterSize3Photon, int(ref_x):int(ref_x)+clusterSize3Photon]
_hists['h1_3PhotonClusterEnergy'].Fill(np.sum(cluster_3Photon))
for i in range(clusterSize3Photon):
for j in range(clusterSize3Photon):
_hists['h2_3PhotonSumSample'].Fill(i, j, cluster_3Photon[i, j])
if 'writeClusters' in _cfg and _cfg['writeClusters'] == True:
cluster_3Photon_list.append(cluster_3Photon)
refpoint_3Photon_list.append([int(ref_x)+Roi[0], int(ref_y)+Roi[2]])
_hists['h2_ClusterSizeVsEnergy'].Fill(energy, clusterSize)
_diameterX = max(xs) - min(xs) + 1
_diameterY = max(ys) - min(ys) + 1