From 6735f41390fe58dd2a50d0a79a884230b97efd37 Mon Sep 17 00:00:00 2001 From: Anna Bergamaschi Date: Thu, 21 Mar 2024 16:36:38 +0100 Subject: [PATCH] added some python script examples --- examples/cluster_example.py | 3 +- examples/clustersFunctions.py | 15 ++++++ examples/color_imaging.py | 47 ++++++++++++++++++ examples/readClusters.py | 89 +++++++++++++++++++++++++++++++++++ 4 files changed, 153 insertions(+), 1 deletion(-) create mode 100644 examples/color_imaging.py create mode 100644 examples/readClusters.py diff --git a/examples/cluster_example.py b/examples/cluster_example.py index 1fd4749..4593434 100644 --- a/examples/cluster_example.py +++ b/examples/cluster_example.py @@ -17,7 +17,7 @@ emax=30 ecutmin=8 ecutmax=12 etabins=251 -csize=2 +csize=3 gain=150 nbins=100 indmin=1 @@ -78,6 +78,7 @@ for i in range(indmin,indmax): fig, ax = plt.subplots() ax.plot(ebins[:-1],sp) +ax.plot(ebins[:-1],spff) #ax.set_yscale('log') fig.show() """ diff --git a/examples/clustersFunctions.py b/examples/clustersFunctions.py index b08ae4b..8709318 100644 --- a/examples/clustersFunctions.py +++ b/examples/clustersFunctions.py @@ -49,6 +49,21 @@ for i in range(1,21): print(ff) r = cr.ClusterFileReader(ff) """ +def color_images(r, emin, emax, ebins, xmin, xmax, ymin, ymax, images=None, csize=3,gain=150): + n=0 + while (cl:=r.read(100000,None)).size: + v=cr.clusterize(csize,cl['data']) + vv=[cl['x'],cl['y'],v['tot']/gain] + image,bins=np.histogramdd(vv, bins=[xmax-xmin,yxmax-ymin,ebins], range=[[xmin,xmax-1],[ymin,ymax-1],[emin,emax]]) + + if images is None: + images = image.copy() + else: + images=images+image + n+=v.shape[0] + #print(v.shape) + print("Read ",n," clusters",np.sum(images)) + return images, bins def analyze_clusters(r, emin, emax, ecutmin, ecutmax, xmin, xmax, ymin, ymax, ietax=None, ietay=None, im=None, sp=None, etas=None, intim=None, csize=3,gain=150, nbins=100, etabins=250, subpix=5): diff --git a/examples/color_imaging.py b/examples/color_imaging.py new file mode 100644 index 0000000..5ad9f9d --- /dev/null +++ b/examples/color_imaging.py @@ -0,0 +1,47 @@ +import os, sys +from pathlib import Path +import matplotlib.pyplot as plt +import numpy as np +#from creader import ClusterFileReader +import creader as cr +import clustersFunctions as cf + +fname = "/mnt/myData/230914_30s_star_100um_nofi/star_" +xmin=0 +xmax=400 +ymin=0 +ymax=400 +emin=0 +emax=30 + +csize=3 +gain=150 +nbins=100 +indmin=1 +indmax=20 + + + +im=None + + + +for i in range(indmin,indmax+1): + #ff=fname + ff=fname+str(i)+".clust" + print(ff) + r = cr.ClusterFileReader(ff) + im, bins=cf.color_images(r, emin, emax, nbins, xmin, xmax, ymin, ymax, im, csize, gain) + +sp=np.sum(im,axis=(0,1)) + +fig, ax = plt.subplots() +ax.plot(bins[2][:-1],sp) +##ax.set_yscale('log') +fig.show() + +cu=np.sum(im[:,:,10:35],axis=2) +cf.plot_colz(cu) + +mo=np.sum(im[:,:,50:70],axis=2) +cf.plot_colz(mo) diff --git a/examples/readClusters.py b/examples/readClusters.py new file mode 100644 index 0000000..cfbe6b7 --- /dev/null +++ b/examples/readClusters.py @@ -0,0 +1,89 @@ +import numpy as np +from numpy.lib import recfunctions as rfn +import matplotlib.pyplot as plt +import sys + +#energyStep=100 +energyMax=40000 +energyBins=400 +clusterSize=3 + +dtypeCluster = [('frameNr', np.int32),('coord', (np.int16,2)),('data',(np.int32,clusterSize*clusterSize))] + +def read_cluster_file(filename): + fd = open(filename,'rb') + clusters = np.fromfile(fd,dtype=dtypeCluster) + fd.close + #rfn.drop_fields(clusters, 'frameNr') + uCl = rfn.structured_to_unstructured(clusters) + return uCl + +def getEnergyArray(data): + off=3 + enArray = np.sum(data[:,off:], axis=-1) + print(data.shape) + Q = np.empty((4,data.shape[0])) + print(Q.shape, data.shape) + Q[0,:]=data[:,off+0]+data[:,off+1]+data[:,off+3]+data[:,off+4] + Q[1,:]=data[:,off+1]+data[:,off+2]+data[:,off+4]+data[:,off+5] + Q[2,:]=data[:,off+3]+data[:,off+4]+data[:,off+6]+data[:,off+7] + Q[3,:]=data[:,off+4]+data[:,off+5]+data[:,off+7]+data[:,off+8] + print(Q) + quadArray = np.max(Q,axis=0) + return enArray,quadArray + +def spectrum_roi(x,y,en,xmin,xmax,ymin,ymax): + global energyBins + global energymax + roi=np.where(np.logical_and(np.logical_and(np.logical_and(x>=xmin,x<=xmax),y>=ymin),y<=ymax)) + spectrum3,xedges2=np.histogram(en[roi],bins=energyBins,range=(0,energyMax)) + return spectrum3,xedges2[1:] + +def image_cut(x,y,en,emin,emax): + energyCut=np.where(np.logical_and(en>=emin,en<=emax)) + print(x[energyCut].shape,y[energyCut].shape) + image,xedges,yedges=np.histogram2d(x[energyCut],y[energyCut],bins=400) + return image + + +fname="/mnt/moench_data/tomcat_fluorescence_24022020/clusters_tomo/cu_fibers_27keV_17.clust" +if len(sys.argv)>1: + fname=sys.argv[1] + +cl=read_cluster_file(fname) +print("file read") + +enADC,quadADC=getEnergyArray(cl) +print("energy array done",enADC,quadADC) + +en=enADC*1000./150. +quad=quadADC*1000./150. +print("energy conversion done:",en,quad) + +x=cl[:,1] +y=cl[:,2] + +#him=get_hyperimage(x,y,cl) +#print("hyperimage done") +xmin=0 +xmax=400 +ymin=00 +ymax=400 + + +spectrum2,xedges1= spectrum_roi(x,y,quad,xmin,xmax,ymin,ymax) +fig2, axs2 = plt.subplots() +axs2.plot(xedges1,spectrum2,"b-") +fig2.show() +print("sp plotted") + +emin=0 +emax=40000 +image=image_cut(x,y,quad,emin,emax) +fig, axs = plt.subplots() +v=axs.imshow(image,vmax=np.mean(image)*5*np.sqrt(np.var(image)),origin='upper',cmap=plt.cm.jet) +fig.colorbar(v, ax=axs) +fig.show() + +print("To get the spectrum for a certain ROI use: sp,x=spectrum_roi(x,y,en,xmin,xmax,ymin,ymax)") +print("To the get the image with a certain energy cut use im=image_cut(x,y,en,emin,emax)")