Untitled

mail@pastecode.io avatar
unknown
plain_text
a month ago
19 kB
2
Indexable
Never
from PyQt5 import QtWidgets, uic, QtCore
import sys
from superqt import QLabeledRangeSlider

#Install ImageJ-PlugIn: EPICS AreaDetector NTNDA-Viewer, look for the channel specified here under channel_name, consider multiple users on servers!!!
import numpy as np

channel_name = 'BAMline:Radio_Theta2_avg'
standard_path = "/raid/CT/2023/2023_03/Markoetter/daisy_lam3_top_mystery/" # '/mnt/raid/CT/2022/'

import numpy
from PIL import Image
import h5py
import tomopy
import math
import time
import os
import cv2
import scipy
import pvaccess as pva      #to install search for "pvapy"
from signalalignment import *
import time
from superqt import QRangeSlider

position = (0, 1000)
fontScale = 5
fontColor = (0, 0, 0)
lineType = 1
thickness = 10

# Function to crop arrays to the same size along the z-dimension
def crop_arrays_to_same_size(arr1, arr2, signal):
    z1, y1, x1 = arr1.shape
    z2, y2, x2 = arr2.shape

    # Check if the arrays have the same size in the y and x dimensions
    if y1 != y2 or x1 != x2:
        raise ValueError("The arrays must have the same size in the y and x dimensions.")

    # Crop the arrays along the z-dimension if needed, from the beginning of the array if signal is true, from the end if signal is false
    if signal:
        if z1 > z2:
            print('cropping projections from the beginning')
            arr1 = arr1[(z1 - z2):, :, :]
        elif z2 > z1:
            print('cropping flat fields from the beginning')
            arr2 = arr2[(z2 - z1):, :, :]
        else:
            print('no cropping needed')
    else:
        if z1 > z2:
            print('cropping projections from the end')
            arr1 = arr1[:(z2 - z1), :, :]
        elif z2 > z1:
            print('cropping flat fields from the end')
            arr2 = arr2[:(z1 - z2), :, :]
        else:
            print('no cropping needed')

    return arr1, arr2

class Radio_Theta2_avg(QtWidgets.QMainWindow):
    def __init__(self):
        super(Radio_Theta2_avg, self).__init__()
        uic.loadUi('Radio_Theta2_avg_test.ui', self)
        self.setWindowTitle('Radio_Theta2_avg')

        self.FFslider = QLabeledRangeSlider(QtCore.Qt.Orientation.Horizontal)
        self.projslider = QLabeledRangeSlider(QtCore.Qt.Orientation.Horizontal)
        #connect buttons to actions
        self.pushLoad.clicked.connect(self.set_path)
        self.pushCompute.clicked.connect(self.compute)
        self.doubleSpinBox_Shift.valueChanged.connect(self.check)
        self.x_start.valueChanged.connect(self.crop_changed)
        self.x_end.valueChanged.connect(self.crop_changed)
        self.y_start.valueChanged.connect(self.crop_changed)
        self.y_end.valueChanged.connect(self.crop_changed)
        self.z_start.valueChanged.connect(self.crop_changed)
        self.z_end.valueChanged.connect(self.crop_changed)
        self.projslider.valueChanged.connect(self.crop_changed)
        self.FFslider.valueChanged.connect(self.crop_changed)

        self.doubleSpinBox_Shift.valueChanged.connect(self.reset_table)

        #### from tomostream.py, nikitinvv git, micha
        # pva type channel that contains projection and metadata
        #self.pva_structure = pva.Channel('PCO1600:Pva1:Image')

        # create pva type pv for reconstruction by copying metadata from the data pv, but replacing the sizes
        # This way the ADViewer (NDViewer) plugin can be also used for visualizing reconstructions.
        #pva_image_data = self.pva_structure.get('')
        #pva_image_dict = pva_image_data.getStructureDict()
        pva_image_dict = {'value': ({'booleanValue': [pva.pvaccess.ScalarType.BOOLEAN], 'byteValue':
            [pva.pvaccess.ScalarType.BYTE], 'shortValue': [pva.pvaccess.ScalarType.SHORT], 'intValue':
            [pva.pvaccess.ScalarType.INT], 'longValue': [pva.pvaccess.ScalarType.LONG], 'ubyteValue':
            [pva.pvaccess.ScalarType.UBYTE], 'ushortValue': [pva.pvaccess.ScalarType.USHORT], 'uintValue':
            [pva.pvaccess.ScalarType.UINT], 'ulongValue': [pva.pvaccess.ScalarType.ULONG], 'floatValue':
            [pva.pvaccess.ScalarType.FLOAT], 'doubleValue': [pva.pvaccess.ScalarType.DOUBLE]},), 'codec':
            {'name': pva.pvaccess.ScalarType.STRING, 'parameters': ()}, 'compressedSize':
            pva.pvaccess.ScalarType.LONG, 'uncompressedSize': pva.pvaccess.ScalarType.LONG, 'dimension':
            [{'size': pva.pvaccess.ScalarType.INT, 'offset': pva.pvaccess.ScalarType.INT, 'fullSize':
                pva.pvaccess.ScalarType.INT, 'binning': pva.pvaccess.ScalarType.INT, 'reverse':
                pva.pvaccess.ScalarType.BOOLEAN}], 'uniqueId': pva.pvaccess.ScalarType.INT, 'dataTimeStamp':
            {'secondsPastEpoch': pva.pvaccess.ScalarType.LONG, 'nanoseconds': pva.pvaccess.ScalarType.INT,
             'userTag': pva.pvaccess.ScalarType.INT}, 'attribute':
            [{'name': pva.pvaccess.ScalarType.STRING, 'value': (), 'descriptor': pva.pvaccess.ScalarType.STRING,
              'sourceType': pva.pvaccess.ScalarType.INT, 'source': pva.pvaccess.ScalarType.STRING}], 'descriptor':
            pva.pvaccess.ScalarType.STRING, 'alarm': {'severity': pva.pvaccess.ScalarType.INT, 'status':
            pva.pvaccess.ScalarType.INT, 'message': pva.pvaccess.ScalarType.STRING}, 'timeStamp':
            {'secondsPastEpoch': pva.pvaccess.ScalarType.LONG, 'nanoseconds': pva.pvaccess.ScalarType.INT, 'userTag':
                pva.pvaccess.ScalarType.INT}, 'display': {'limitLow': pva.pvaccess.ScalarType.DOUBLE, 'limitHigh':
            pva.pvaccess.ScalarType.DOUBLE, 'description': pva.pvaccess.ScalarType.STRING, 'format':
            pva.pvaccess.ScalarType.STRING, 'units': pva.pvaccess.ScalarType.STRING}}

        self.pv_rec = pva.PvObject(pva_image_dict)
        self.pvaServer = pva.PvaServer(channel_name, self.pv_rec)
        self.Qchannel_name.setText(channel_name)
        self.pvaServer.start()
        self.normeverythingtiff.clicked.connect(self.recoeveverythingtiff)
        self.crop_changed = False
        self.shift_calculated = False
        self.crop_beginning_checkbox.setChecked(True)


    def set_path(self):

        self.new = 1
        # ask for hdf5-file #proj
        path_klick_proj = QtWidgets.QFileDialog.getOpenFileName(self,'Select projections hdf5-file, please.', standard_path)
        self.path_klick_proj = path_klick_proj[0]
        print('path klicked: ', self.path_klick_proj)

        # ask for hdf5-file #ff
        path_klick_ff = QtWidgets.QFileDialog.getOpenFileName(self, 'Select flat-field hdf5-file, please.', standard_path)
        self.path_klick_ff = path_klick_ff[0]
        print('path klicked: ', self.path_klick_ff)

        self.Qfilename.setText(self.path_klick_proj)

        # link a volume to the proj hdf-file
        f_proj = h5py.File(self.path_klick_proj, 'r')
        self.vol_proxy_proj = f_proj['/entry/data/data']
        print('raw data proj volume size: ', self.vol_proxy_proj.shape)

        # link a volume to the ff hdf-file
        f_ff = h5py.File(self.path_klick_ff, 'r')
        self.vol_proxy_ff = f_ff['/entry/data/data']
        #self.vol_proxy_ff = np.divide(self.vol_proxy_ff,self.vol_proxy_ff)
        print('raw data ff volume size: ', self.vol_proxy_ff.shape)

        self.x_start.setMinimum(0)
        self.x_start.setMaximum(self.vol_proxy_proj.shape[1])
        self.x_end.setMaximum(self.vol_proxy_proj.shape[1])
        self.x_end.setValue(self.vol_proxy_proj.shape[1])

        self.y_start.setMinimum(0)
        self.y_start.setMaximum(self.vol_proxy_proj.shape[2])
        self.y_end.setMaximum(self.vol_proxy_proj.shape[2])
        self.y_end.setValue(self.vol_proxy_proj.shape[2])

        self.height_alignment_slice.setMinimum(0)
        self.height_alignment_slice.setMaximum(self.vol_proxy_proj.shape[2])
        self.height_alignment_slice.setValue(round(self.vol_proxy_proj.shape[2]/2))
        self.norm_then_divide_checkbox.setChecked(False)

        self.Proj,self.FFs = crop_arrays_to_same_size(self.vol_proxy_proj,self.vol_proxy_ff,self.crop_beginning_checkbox.isChecked())

        self.min_size = min(self.FFs.shape[0],self.Proj.shape[0])

        self.FFslider.setRange(0, self.vol_proxy_ff.shape[0])
        self.projslider.setRange(0, self.vol_proxy_proj.shape[0])



        self.sliders_layout.layout().addWidget(QtWidgets.QLabel('FF'))
        self.sliders_layout.layout().addWidget(self.FFslider)
        self.sliders_layout.layout().addWidget(QtWidgets.QLabel('Projections'))
        self.sliders_layout.layout().addWidget(self.projslider)

        self.z_start.setMaximum(self.min_size)
        self.z_end.setMaximum(self.min_size)
        self.z_end.setValue(self.min_size)
        self.z_end.setMinimum(1)
        self.z_start.setMinimum(0)
        self.z_start.setValue(0)

        self.load()

    def set_path2(self):
        self.even_files = []
        self.odd_files = []

        path_klick_proj = QtWidgets.QFileDialog.getOpenFileName(self, 'Select folder containing all sub-folders please', standard_path)
        self.path_klick_proj = path_klick_proj[0]
        self.folder_names = [f for f in os.listdir(self.path_klick_proj) if os.path.isdir(os.path.join(self.path_klick_proj, f))]
        self.sorted_folder_names = sorted(self.folder_names, key=lambda x: int(x.split("_")[0]))

        for folder in self.sorted_folder_names:
            number = folder.split("_")[0]
    def crop_changed(self):
        self.crop_changed = True

    #function to change values in the second row of the table. the subsequent values after the one changed also change, with a delta of the value changed
    def table_changed(self, item):
        if item.row() == 1:
            old_value = self.saved_table[item.column()]
            delta = float(item.text()) - old_value
            self.tablewidget_shifts.disconnect()
            for i in range(item.column(), self.tablewidget_shifts.columnCount()):
                self.tablewidget_shifts.item(1,i).setText(str(self.saved_table[i] + delta))
        self.save_table()
        self.tablewidget_shifts.itemChanged.connect(self.table_changed)

    def save_table(self):
        self.saved_table = []  # This list will hold the values of the second row

        for column in range(self.tablewidget_shifts.columnCount()):
            item = self.tablewidget_shifts.item(1, column)  # Row index 1 since it's the second row (0-indexed)
            if item is not None:  # Check if the item exists
                self.saved_table.append(float(item.text()))

    def load(self):
        print('Projection taken from ', self.projslider.value()[0], ' to ', self.projslider.value()[1])
        print('FF taken from ', self.FFslider.value()[0], ' to ', self.FFslider.value()[1])

        self.loadedProj, self.loadedFFs = crop_arrays_to_same_size(self.vol_proxy_proj[
                                                       self.projslider.value()[0]:self.projslider.value()[1],
                                                       self.x_start.value():self.x_end.value(),
                                                       self.y_start.value():self.y_end.value()]-self.spinBox_DF.value(),
                                                        self.vol_proxy_ff[
                                                        self.FFslider.value()[0]:self.FFslider.value()[1],
                                                        self.x_start.value():self.x_end.value(),
                                                        self.y_start.value():self.y_end.value()]-self.spinBox_DF.value(),
                                                       self.crop_beginning_checkbox.isChecked())


        self.shifts = numpy.zeros((self.FFs.shape[0],2))

        print('proj shape: ', self.loadedProj.shape)
        print('ff shape: ', self.loadedFFs.shape)
        print('shifts: ', self.shifts.shape)

        self.Norm_stack = numpy.zeros((self.loadedProj.shape[0],self.loadedProj.shape[1],self.loadedProj.shape[2]))
        self.pv_rec['dimension'] = [
            {'size': self.Norm_stack.shape[2], 'fullSize': self.Norm_stack.shape[2], 'binning': 1},
            {'size': self.Norm_stack.shape[1], 'fullSize': self.Norm_stack.shape[1], 'binning': 1}]

        print('volume loaded: ', self.Norm_stack.shape)
        #feed the line edit with the the manual shift array
        self.tablewidget_shifts.setRowCount(2)
        self.tablewidget_shifts.setColumnCount(self.Norm_stack.shape[0])

        #fill the manual shift array with the pixel/image shift values, i.e pxl/image * img number
        self.manual_shift = numpy.arange(self.Norm_stack.shape[0])*self.doubleSpinBox_Shift.value()

        #we fill the table widget with the manual shift array, 2 rows, one for the index, one for the shift, and norm.stack.shape[0] columns
        for i in range(self.Norm_stack.shape[0]):
            print(i, str(self.manual_shift[i]))
            self.tablewidget_shifts.setItem(0,i,QtWidgets.QTableWidgetItem(str(i)))
            self.tablewidget_shifts.setItem(1,i,QtWidgets.QTableWidgetItem(str(self.manual_shift[i])))

        self.save_table()
        self.stop_changing_table_signal = True
        self.tablewidget_shifts.itemChanged.connect(self.table_changed)

    def reset_table(self):
        self.manual_shift = numpy.arange(self.Norm_stack.shape[0]) * self.doubleSpinBox_Shift.value()
        for i in range(self.Norm_stack.shape[0]):
            self.tablewidget_shifts.setItem(0, i, QtWidgets.QTableWidgetItem(str(i)))
            self.tablewidget_shifts.setItem(1, i, QtWidgets.QTableWidgetItem(str(self.manual_shift[i])))

        self.save_table()
        self.stop_changing_table_signal = True
        self.tablewidget_shifts.itemChanged.connect(self.table_changed)

    def check(self):
        if self.checkBox_auto.isChecked():
           self.compute()

    def compute(self):
        self.pushCompute.setText('Busy')
        self.pushLoad.setEnabled(False)
        self.pushCompute.setEnabled(False)  # DOES NOT SEEM TO WORK!

        if self.crop_changed:
            self.load()
            self.crop_changed = False

        self.Norm_stack2 = self.Norm_stack.copy()
        self.Proj_shift = numpy.zeros((self.Norm_stack.shape[0],self.Norm_stack.shape[1],self.Norm_stack.shape[2]))
        self.FF_shift = numpy.zeros((self.Norm_stack.shape[0],self.Norm_stack.shape[1],self.Norm_stack.shape[2]))
        self.number_slices = self.Norm_stack.shape[0]

        i = 0

        while i < self.number_slices:
            print(i, ' of ',self.number_slices)
            if self.norm_then_divide_checkbox.isChecked():
                print('Averaging then dividing')
                self.proj_shift = float(self.tablewidget_shifts.item(1,i).text())
                print('Image number: ', i, ' shift: ', self.proj_shift)
                self.Proj_shift[i,:,:] = scipy.ndimage.shift(self.loadedProj[i,:,:],
                                                              (self.proj_shift,0),
                                                              order=int(self.intSpinbox_Interpolationorder.value()),
                                                              mode=self.interpolation_mode_combobox.currentText(),
                                                              prefilter=True)
                self.FF_shift[i,:,:] = scipy.ndimage.shift(self.loadedFFs[i,:,:],
                                                                (self.proj_shift,0),
                                                                order=int(self.intSpinbox_Interpolationorder.value()),
                                                                mode=self.interpolation_mode_combobox.currentText(),
                                                                prefilter=True)
            else:
                print('Dividing then averaging')
                self.proj_shift = float(self.tablewidget_shifts.item(1,i).text())
                #let's try recalculating the norm stack and dividing the shifted loadedProj by the shifted loadedFFs
                print('Image number: ', i, ' shift: ', self.proj_shift)
                self.Norm_stack2[i,:,:] = np.divide(scipy.ndimage.shift(self.loadedProj[i,:,:],
                                                              (self.proj_shift,0),
                                                              order=int(self.intSpinbox_Interpolationorder.value()),
                                                              mode=self.interpolation_mode_combobox.currentText(),
                                                              prefilter=True),
                                                    scipy.ndimage.shift(self.loadedFFs[i,:,:],
                                                                (self.proj_shift,0),
                                                                order=int(self.intSpinbox_Interpolationorder.value()),
                                                                mode=self.interpolation_mode_combobox.currentText(),
                                                                prefilter=True))
                if self.checkBox_intermediate_images.isChecked():
                    self.pv_rec['value'] = ({'floatValue': np.asarray(
                        self.Norm_stack2[i,:,:], dtype=np.float32).flatten()},)
                    time.sleep(0.1)

            i=i+1

        if self.norm_then_divide_checkbox.isChecked():
            self.avg = np.divide(np.mean(self.Proj_shift, axis=0),np.mean(self.FF_shift, axis=0), dtype=numpy.float32)
        else:
            self.avg = np.mean(self.Norm_stack2, axis=0, dtype=numpy.float32)

        #write result to pv
        self.pv_rec['value'] = ({'floatValue': self.avg.flatten()},)

        self.meanhigh = numpy.mean(self.avg[717:817,848:848+139])
        self.meanlow = numpy.mean(self.avg[717:817,618:618+139])
        self.stdhigh = numpy.std(self.avg[717:817,848:848+139])
        self.stdlow = numpy.std(self.avg[717:817,618:618+139])
        self.snrtext.setText(str((abs(self.meanhigh-self.meanlow))/(self.stdhigh**2+self.stdlow**2)**0.5))

        self.pushCompute.setText('Compute')
        self.pushLoad.setEnabled(True)
        self.pushCompute.setEnabled(True)

    def recoeveverythingtiff(self):
        self.pushCompute.setText('Busy')
        self.pushLoad.setEnabled(False)
        self.pushCompute.setEnabled(False)

# #no idea why we need this, but it wouldn't work without it ;-)
if __name__ == "__main__":
    import sys
    app = QtWidgets.QApplication(sys.argv)
    main = Radio_Theta2_avg()
    main.show()
    sys.exit(app.exec_())
Leave a Comment