From 19575b6b8b7d3a34fde89df265bb77639cc17a98 Mon Sep 17 00:00:00 2001 From: nducros Date: Fri, 27 Jan 2023 19:42:38 +0100 Subject: [PATCH 01/11] Update reconstruct & add reorder and subsample --- scripts/reconstruction_neural_network.py | 66 +++++ spas/__init__.py | 4 +- spas/reconstruction_nn.py | 359 +++++++++++++++-------- spas/visualization.py | 7 +- 4 files changed, 306 insertions(+), 130 deletions(-) create mode 100644 scripts/reconstruction_neural_network.py diff --git a/scripts/reconstruction_neural_network.py b/scripts/reconstruction_neural_network.py new file mode 100644 index 0000000..9c1b786 --- /dev/null +++ b/scripts/reconstruction_neural_network.py @@ -0,0 +1,66 @@ +# -*- coding: utf-8 -*- +""" +Created on Thu Jan 26 16:49:14 2023 + +@author: ducros +""" + +#%% +from spas import ReconstructionParameters, setup_reconstruction + +network_param = ReconstructionParameters( + # Reconstruction network + M = 64*64, # Number of measurements + img_size = 128, # Image size + arch = 'dc-net', # Main architecture + denoi = 'unet', # Image domain denoiser + subs = 'rect', # Subsampling scheme + + # Training + data = 'imagenet', # Training database + N0 = 10, # Intensity (max of ph./pixel) + + # Optimisation (from train2.py) + num_epochs = 30, # Number of training epochs + learning_rate = 0.001, # Learning Rate + step_size = 10, # Scheduler Step Size + gamma = 0.5, # Scheduler Decrease Rate + batch_size = 256, # Size of the training batch + regularization = 1e-7 # Regularisation Parameter + ) + +cov_path = '../../stat/ILSVRC2012_v10102019/Cov_8_128x128.npy' +model_folder = '../../model_v2/' + +model, device = setup_reconstruction(cov_path, model_folder, network_param) + +#%% Load data and meta data +import numpy as np + +# data +data_path = '../../pilot-spihim/cat_linear/Cat_LinearColoredFilter' +full_path = data_path + '_spectraldata.npz' +meas = np.load(full_path)['spectral_data'] + +#%% Spectral binning +from spas import read_metadata, spectral_binning + +# load meta data to get wavelength +meta_path = data_path + '_metadata.json' +_, acquisition_param, _, _ = read_metadata(meta_path) +wavelengths = acquisition_param.wavelengths + +# bin raw data +meas_bin, wavelengths_bin, _ = spectral_binning(meas.T, wavelengths, 530, 730, 4) + +#%% Reorder and subsample +from spas.reconstruction_nn import reorder_subsample +meas_bin_2 = reorder_subsample(meas_bin, acquisition_param, network_param) + +#%% Reconstruct +from spas import reconstruct +rec = reconstruct(model, device, meas_bin_2, 1) + +#%% Plot +from spas import plot_color +plot_color(rec, wavelengths_bin) \ No newline at end of file diff --git a/spas/__init__.py b/spas/__init__.py index f83ed88..a711ef3 100644 --- a/spas/__init__.py +++ b/spas/__init__.py @@ -1,11 +1,11 @@ # -*- coding: utf-8 -*- __author__ = 'Guilherme Beneti Martins' -from .acquisition import init, setup, acquire, disconnect #, init_2arms, disconnect_2arms +#from .acquisition import init, setup, acquire, disconnect #, init_2arms, disconnect_2arms from .metadata import MetaData, AcquisitionParameters from .metadata import DMDParameters, SpectrometerParameters from .metadata import read_metadata, save_metadata -from .generate import * +#from .generate import * from .reconstruction import * from .visualization import * from .noise import * diff --git a/spas/reconstruction_nn.py b/spas/reconstruction_nn.py index d8ea798..81a3999 100644 --- a/spas/reconstruction_nn.py +++ b/spas/reconstruction_nn.py @@ -15,13 +15,28 @@ from multiprocessing import Queue from pathlib import Path +import math import torch import numpy as np from matplotlib import pyplot as plt -from spyrit.learning.model_Had_DCAN import compNet, noiCompNet, DenoiCompNet -from spyrit.learning.nets import load_net + +from spyrit.core.reconstruction import Pinv_Net, DC2_Net + +from spyrit.core.training import load_net +from spyrit.misc.statistics import Cov2Var +from spyrit.misc.walsh_hadamard import walsh2_matrix + +from spyrit.core.Acquisition import Acquisition_Poisson_approx_Gauss +from spyrit.core.Forward_Operator import Forward_operator_Split_ft_had +from spyrit.core.Preprocess import Preprocess_Split_diag_poisson +from spyrit.core.Data_Consistency import Generalized_Orthogonal_Tikhonov, Pinv_orthogonal +from spyrit.core.neural_network import Unet, Identity +from spyrit.misc.sampling import Permutation_Matrix + + from spas.noise import noiseClass +from spas.metadata import AcquisitionParameters class netType(IntEnum): """Possible model architectures. @@ -36,135 +51,216 @@ class netType(IntEnum): class ReconstructionParameters: """Reconstruction parameters for loading a reconstruction model. """ - img_size: int - CR: int - denoise: bool - epochs: int - learning_rate: float - step_size: int - gamma: float - batch_size: int - regularization: float - N0: float - sig: float - - arch_name: InitVar[str] - - _net_arch: int = field(init=False) - - - def __post_init__(self, arch_name): - self.arch_name = arch_name - - - @property - def arch_name(self): - return netType(self._net_arch).name - - - @arch_name.setter - def arch_name(self, arch_name): - self._net_arch = int(netType[arch_name]) + # Reconstruction network + M: int # Number of measurements + img_size: int # Image size + arch: str # Main architecture + denoi: str # Image domain denoiser + subs: str # Subsampling scheme + + # Training + data: str # Training database + N0: float # Intensity (max of ph./pixel) + #sig: float + + # Optimisation (from train2.py) + num_epochs: int # Number of training Epochs + learning_rate: float # Learning Rate + step_size: int # Scheduler Step Size + gamma: float # Scheduler Decrease Rate + batch_size: int # Size of the training batch + regularization: float # Regularisation Parameter + #checkpoint_model: str # Optional path to checkpoint model + #checkpoint_interval: int# Interval between saving model checkpoints -def setup_reconstruction(cov_path: str, mean_path: str, H: np.ndarray, - model_root: str, network_params: ReconstructionParameters - ) -> Tuple[Union[compNet, noiCompNet, DenoiCompNet], str]: +def setup_reconstruction(cov_path: str, + model_folder: str, + network_params: ReconstructionParameters + ) -> Tuple[Union[Pinv_Net,DC2_Net], str]: """Loads a neural network for reconstruction. + Limited to measurements from patterns of size 2**K for reconstruction at + size 2**L, with L > K (e.g., measurements at size 64, reconstructions at + size 128). + Args: cov_path (str): - Path to the covariance matrix. - mean_path (str): - Path to the mean matrix. - H (nd.array): - Hadamard matrix with patterns. - model_root (str): + Path to the covariance matrix used for reconstruction. + model_folder (str): Folder containing trained models for reconstruction. network_params (ReconstructionParameters): Parameters used to load the model. Returns: - Tuple[Union[compNet, noiCompNet, DenoiCompNet], str]: - model (compNet, noiCompNet, DenoiCompNet): + Tuple[Union[Pinv_Net, DC2_Net], str]: + model (Pinv_Net, DC2_Net): Loaded model. device (str): Device to which the model was loaded. """ - net_type = ['c0mp', 'comp', 'pinv', 'free'] - device = torch.device("cuda:0" if torch.cuda.is_available() else "cpu") print(f'Device: {device}') - Cov_had = np.load(cov_path) / network_params.img_size**2 - Mean_had = np.load(mean_path) / network_params.img_size - - suffix = '_N_{}_M_{}_epo_{}_lr_{}_sss_{}_sdr_{}_bs_{}_reg_{}'.format( - network_params.img_size, network_params.CR, - network_params.epochs, network_params.learning_rate, - network_params.step_size, network_params.gamma, - network_params.batch_size, network_params.regularization) - - recon_type = "" - if network_params.N0 == 0: - train_type = '' - else: - train_type = '_N0_{}_sig_{}'.format(network_params.N0, - network_params.sig) - if network_params.denoise == True: - recon_type+="_Denoi" - - # Training parameters - arch = network_params.arch_name - - suffix = 'NET_' + arch + train_type + recon_type + suffix - - title = Path(model_root) / suffix - - if network_params.N0 == 0: - model = compNet( - network_params.img_size, - network_params.CR, - Mean_had, - Cov_had, - network_params._net_arch, - H) + Cov_rec = np.load(cov_path) + H = walsh2_matrix(network_params.img_size) - elif network_params.denoise == False: - model = noiCompNet( - network_params.img_size, - network_params.CR, - Mean_had, - Cov_had, - network_params._net_arch, - network_params.N0, - network_params.sig, - H) - - elif network_params.denoise == True: - model = DenoiCompNet( - network_params.img_size, - network_params.CR, - Mean_had, - Cov_had, - network_params._net_arch, - network_params.N0, - network_params.sig, - H, - None) - - torch.cuda.empty_cache() - + # Rectangular sampling + # N.B.: Only for measurements from patterns of size 2**K reconstructed at + # size 2**L, with L > K (e.g., measurements are at size 64, reconstructions + # at size 128. + Ord = np.ones((network_params.img_size, network_params.img_size)) + n_sub = math.ceil(network_params.M**0.5) + Ord[:,n_sub:] = 0 + Ord[n_sub:,:] = 0 + + # Init network + Perm_rec = Permutation_Matrix(Ord) + Hperm = Perm_rec @ H + Pmat = Hperm[:network_params.M,:] + + # init + Forward = Forward_operator_Split_ft_had(Pmat, Perm_rec, + network_params.img_size, + network_params.img_size) + Noise = Acquisition_Poisson_approx_Gauss(network_params.N0, Forward) + Prep = Preprocess_Split_diag_poisson(network_params.N0, + network_params.M, + network_params.img_size**2) + + Denoi = Unet() + Cov_perm = Perm_rec @ Cov_rec @ Perm_rec.T + DC = Generalized_Orthogonal_Tikhonov(sigma_prior = Cov_perm, + M = network_params.M, + N = network_params.img_size**2) + model = DC2_Net(Noise, Prep, DC, Denoi) + + # load + net_folder = '{}_{}_{}'.format( + network_params.arch, network_params.denoi, + network_params.data) + + suffix = '_{}_N0_{}_N_{}_M_{}_epo_{}_lr_{}_sss_{}_sdr_{}_bs_{}_reg_{}'.format( + network_params.subs, network_params.N0, + network_params.img_size, network_params.M, + network_params.num_epochs, network_params.learning_rate, + network_params.step_size, network_params.gamma, + network_params.batch_size, network_params.regularization) + + torch.cuda.empty_cache() # keep this here? + title = Path(model_folder) / net_folder / (net_folder + suffix) load_net(title, model, device) model = model.to(device) return model, device +def reorder_subsample(meas: np.ndarray, + acqui_param: AcquisitionParameters, + recon_param: ReconstructionParameters, + ) -> np.ndarray: + """Reorder and subsample measurements + + Args: + meas (np.ndarray): + Spectral measurements with dimensions (N_wavelength x M_acq), where + M_acq is the number of acquired patterns + acqui_param (AcquisitionParameters): + Parameters used during the acquisition of the spectral measurments + recon_param (ReconstructionParameters): + Parameters of the reconstruction. -def reconstruct(model: Union[compNet, noiCompNet, DenoiCompNet], - device: str, spectral_data: np.ndarray, batches : int, - noise_model: noiseClass = None, is_process: bool = False) -> np.ndarray: + Returns: + (np.ndarray): + Spectral measurements with dimensions (N_wavelength x M_rec), where + M_rec is the number of patterns considered for reconstruction. + Acquisitions can be subsampled a posteriori, leadind to M_rec < M_acq + + """ + # Dimensions (N.B: images are assumed to be square) + N_acq = acqui_param.pattern_dimension_x + N_rec = recon_param.img_size + N_wav = meas.shape[0] + + # Order used for acquisistion + Ord_acq = -np.array(acqui_param.patterns)[::2]//2 # pattern order + Ord_acq = np.reshape(Ord_acq, (N_acq,N_acq)) # sampling map + Perm_acq = Permutation_Matrix(Ord_acq).T + + + # Order used for reconstruction + if recon_param.subs == 'rect': + Ord_rec = np.ones((N_rec, N_rec)) + n_sub = math.ceil(recon_param.M**0.5) + Ord_rec[:,n_sub:] = 0 + Ord_rec[n_sub:,:] = 0 + Perm_rec = Permutation_Matrix(Ord_rec) + + # + meas = meas.T + + # Subsampled acquisition permutation matrix (and zero filling if necessary) + if N_rec > N_acq: + + # Square subsampling in the "natural" order + Ord_sub = np.zeros((N_rec,N_rec)) + Ord_sub[:N_acq,:N_acq]= -np.arange(-N_acq**2,0).reshape(N_acq,N_acq) + Perm_sub = Permutation_Matrix(Ord_sub) + + # Natural order measurements (N_acq resolution) + Perm_raw = np.zeros((2*N_acq**2,2*N_acq**2)) + Perm_raw[::2,::2] = Perm_acq.T + Perm_raw[1::2,1::2] = Perm_acq.T + meas = Perm_raw @ meas + + # Zero filling (needed only when reconstruction resolution is higher + # than acquisition res) + zero_filled = np.zeros((2*N_rec**2, N_wav)) + zero_filled[:2*N_acq**2,:] = meas + + meas = zero_filled + + Perm_raw = np.zeros((2*N_rec**2,2*N_rec**2)) + Perm_raw[::2,::2] = Perm_sub.T + Perm_raw[1::2,1::2] = Perm_sub.T + + meas = Perm_raw @ meas + + elif N_rec == N_acq: + Perm_sub = Perm_acq[:N_rec**2,:].T + + elif N_rec < N_acq: + # Square subsampling in the "natural" order + Ord_sub = np.zeros((N_acq,N_acq)) + Ord_sub[:N_rec,:N_rec]= -np.arange(-N_rec**2,0).reshape(N_rec,N_rec) + Perm_sub = Permutation_Matrix(Ord_sub) + Perm_sub = Perm_sub[:N_rec**2,:] + Perm_sub = Perm_sub @ Perm_acq.T + + #Reorder measurements when reconstruction order is not "natural" + if N_rec <= N_acq: + # Get both positive and negative coefficients permutated + Perm = Perm_rec @ Perm_sub + Perm_raw = np.zeros((2*N_rec**2,2*N_acq**2)) + + elif N_rec > N_acq: + Perm = Perm_rec + Perm_raw = np.zeros((2*N_rec**2,2*N_rec**2)) + + Perm_raw[::2,::2] = Perm + Perm_raw[1::2,1::2] = Perm + meas = Perm_raw @ meas + + return meas[:2*recon_param.M,:].T + + +def reconstruct(model: Union[Pinv_Net, DC2_Net], + device: str, + spectral_data: np.ndarray, + batches : int, + #noise_model: noiseClass = None, + is_process: bool = False + ) -> np.ndarray: """Reconstructs images from experimental data. Using a loaded model, reconstructs images in batches. It can be used in @@ -172,7 +268,7 @@ def reconstruct(model: Union[compNet, noiCompNet, DenoiCompNet], the same time. Args: - model (Union[compNet, noiCompNet, DenoiCompNet]): + model (Union[Pinv_Net,DC2_Net]): Pre-loaded model for reconstruction. device (str): Device to which the model was loaded. Used to make sure the spectral @@ -185,7 +281,7 @@ def reconstruct(model: Union[compNet, noiCompNet, DenoiCompNet], have enough memory to reconstruct all data at a single time. E.g. when reconstructing data distributed along many wavelengths, there may be too many data. - noise_model (noiseClass, optional): + ?? noise_model (noiseClass, optional): Loaded noise model in case the reconstruction should use a denoising method. Defaults to None. is_process (bool, optional): @@ -198,48 +294,57 @@ def reconstruct(model: Union[compNet, noiCompNet, DenoiCompNet], spectral dimension x image size x image size. """ - # Implemented only the case for Denoi reconstruction proportion = spectral_data.shape[0]//batches # Amount of wavelengths per batch - recon = np.zeros((spectral_data.shape[0], 64, 64)) + img_size = model.Acq.FO.h # image assumed to be square + recon = np.zeros((spectral_data.shape[0], img_size, img_size)) start = perf_counter_ns() + model.PreP.set_expe() + with torch.no_grad(): for batch in range(batches): - lambda_indeces = range(proportion * batch, proportion * (batch+1)) + lambda_indices = range(proportion * batch, proportion * (batch+1)) info = (f'batch {batch},' - f'reconstructed wavelength range: {lambda_indeces}') + f'reconstructed wavelength range: {lambda_indices}') if is_process: info = '#Recon process:' + info print(info) - C = noise_model.mu[lambda_indeces] - s = noise_model.sigma[lambda_indeces] - K = noise_model.K[lambda_indeces] + # C = noise_model.mu[lambda_indices] + # s = noise_model.sigma[lambda_indices] + # K = noise_model.K[lambda_indices] - n = len(C) + # n = len(C) - C = torch.from_numpy(C).float().to(device).reshape(n,1,1) - s = torch.from_numpy(s).float().to(device).reshape(n,1,1) - K = torch.from_numpy(K).float().to(device).reshape(n,1,1) + # C = torch.from_numpy(C).float().to(device).reshape(n,1,1) + # s = torch.from_numpy(s).float().to(device).reshape(n,1,1) + # K = torch.from_numpy(K).float().to(device).reshape(n,1,1) - CR = spectral_data.shape[1] + # M = spectral_data.shape[1] + + # torch_img = torch.from_numpy(spectral_data[lambda_indices,:]) + # torch_img = torch_img.float() + # torch_img = torch.reshape(torch_img, (len(lambda_indices), 1, M)) # batches, channels, patterns + # torch_img = torch_img.to(device) + + spectral_data_torch = torch.tensor(spectral_data[lambda_indices,:], + dtype = torch.float, + device = device) - torch_img = torch.from_numpy(spectral_data[lambda_indeces,:]) - torch_img = torch_img.float() - torch_img = torch.reshape(torch_img, (len(lambda_indeces), 1, CR)) # batches, channels, patterns - torch_img = torch_img.to(device) + # spectral_data_torch = \ + # torch.Tensor(spectral_data[lambda_indices,:]).to(device) - torch_recon = model.forward_reconstruct_expe(torch_img, - len(lambda_indeces), 1, model.n, model.n, C, s, K) + recon_torch = model.reconstruct_expe(spectral_data_torch)#, + #len(lambda_indices), 1, model.n, model.n, C, s, K) - recon[lambda_indeces,:,:] = torch_recon.cpu().detach().numpy().squeeze() + recon[lambda_indices,:,:] = recon_torch.cpu().detach().numpy().squeeze() end = perf_counter_ns() @@ -253,7 +358,7 @@ def reconstruct(model: Union[compNet, noiCompNet, DenoiCompNet], return recon -def reconstruct_process(model: Union[compNet, noiCompNet, DenoiCompNet], +def reconstruct_process(model: Union[Pinv_Net, DC2_Net], device: str, queue_to_recon: Queue, queue_reconstructed: Queue, batches: int, noise_model: noiseClass, sleep_time: float = 0.3) -> None: """Performs reconstruction in real-time. diff --git a/spas/visualization.py b/spas/visualization.py index 1c0d082..6e05de7 100644 --- a/spas/visualization.py +++ b/spas/visualization.py @@ -7,7 +7,7 @@ import os from spas import * -from spas.reconstruction_nn import reconstruct +#from spas.reconstruction_nn import reconstruct import numpy as np from matplotlib.colors import ListedColormap from matplotlib import pyplot as plt @@ -405,6 +405,11 @@ def plot_color(F: np.ndarray, wavelengths: np.ndarray, filename: str = None, fig.tight_layout() fig.tight_layout() + + # save + if filename is not None: + fig.savefig(filename) + # plt.show() ######################### IDS CAM visualisationtion ########################### From f18214a5cbecac9c9010a6dcd5c58ba5d3982562 Mon Sep 17 00:00:00 2001 From: nducros Date: Mon, 30 Jan 2023 18:28:55 +0100 Subject: [PATCH 02/11] update doc and switch reconstruction to eval mode --- README.md | 76 ++++++++++++++---------- scripts/reconstruction_neural_network.py | 10 +++- spas/reconstruction_nn.py | 57 ++++++++---------- 3 files changed, 80 insertions(+), 63 deletions(-) diff --git a/README.md b/README.md index a16ac2e..11016f0 100644 --- a/README.md +++ b/README.md @@ -195,7 +195,7 @@ from spas import plot_color plot_color(rec_bin, wavelengths_bin) ``` -* Measurements are saved in the disk (fully sampled) +* Measurements are saved on the disk (fully sampled) Reconstruct the measurements saved as `../meas/my_first_measurement`. @@ -220,7 +220,7 @@ rec = reconstruction_hadamard(acquisition_metadata.patterns, 'walsh', H, meas) ``` ## Reconstruction with a Neural Network -* Measurements are on the disk (fully-sampled here, works too with `pattern_compression=.25`) +* We consider an existing acquisition that was saved on the disk in the `../meas/` folder Read the data: ``` python @@ -235,51 +235,67 @@ from spas import read_metadata, reconstruction_hadamard _, acquisition_parameters, _, _ = read_metadata('../meas/my_first_measurement' + '_metadata.json') ``` -* We have access to a trained network: +* We consider that we have access to a trained network and the covariance matrix associated to it. -An example network can be downloaded [here](https://www.creatis.insa-lyon.fr/~ducros/spyritexamples/2021_ISTE/NET_c0mp_N0_50.0_sig_0.0_Denoi_N_64_M_1024_epo_40_lr_0.001_sss_20_sdr_0.2_bs_256_reg_1e-07.pth) , which you can save in `./models/`. It allows to reconstruction images from only 1024 hadamard coefficients (i.e., 2048 raw measurements): +An example network can be downloaded [here](https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/622b5ea843258e76eab21740). It allows the reconstruction of a 128 x 128 image from only 4096 Hadamard coefficients (i.e., 8192 raw measurements) that correspond to a full acquisition at a 64 x 64 resolution. Its associated covariance matrix can be downloaded [here](https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/63d7f3620386da2747641e1b). ``` python from spas import ReconstructionParameters, setup_reconstruction -network_params = ReconstructionParameters( - img_size=64, - CR=1024, - denoise=True, - epochs=40, - learning_rate=1e-3, - step_size=20, - gamma=0.2, - batch_size=256, - regularization=1e-7, - N0=50.0, - sig=0.0, - arch_name='c0mp',) - -cov_path = '../stats/Cov_64x64.npy' -mean_path = '../stats/Average_64x64.npy' -model_root = '../models/' -import spyrit.misc.walsh_hadamard as wh -H = wh.walsh2_matrix(64)/64 -model, device = setup_reconstruction(cov_path, mean_path, H, model_root, network_params) +network_param = ReconstructionParameters( + # Reconstruction network + M = 64*64, # Number of measurements + img_size = 128, # Image size + arch = 'dc-net', # Main architecture + denoi = 'unet', # Image domain denoiser + subs = 'rect', # Subsampling scheme + + # Training + data = 'imagenet', # Training database + N0 = 10, # Intensity (max of ph./pixel) + + # Optimisation (from train2.py) + num_epochs = 30, # Number of training epochs + learning_rate = 0.001, # Learning Rate + step_size = 10, # Scheduler Step Size + gamma = 0.5, # Scheduler Decrease Rate + batch_size = 256, # Size of the training batch + regularization = 1e-7 # Regularisation Parameter + ) + +cov_path = '../stat/Cov_8_128x128.npy' +model_folder = '../model/' + +model, device = setup_reconstruction(cov_path, model_folder, network_param) ``` -Load noise calibration parameters (provided with the data or computed using tools in `/noise-calibration`) +Load noise calibration parameters (provided with the data or computed using tools in `/noise-calibration`). :warning:: Not used anymore in the current implementation of `spas`. ``` python from spas import load_noise noise = load_noise('../noise-calibration/fit_model.npz') ``` -Bin before reconstruction and plot +Bin the spectral measurements (here, 4 bins between 530 nm and 730 nm) ``` python from spas import spectral_binning -meas_bin, wavelengths_bin, _, noise_bin = spectral_binning(meas.T, acquisition_parameters.wavelengths, 530, 730, 8, noise) +meas_bin, wavelengths_bin, _ = spectral_binning(meas.T, wavelengths, 530, 730, 4) +``` + +Reorder and subsample the spectral measurements +``` python +from spas.reconstruction_nn import reorder_subsample +meas_bin_2 = reorder_subsample(meas_bin, acquisition_param, network_param) +``` + +Reconstruct the spectral images +``` python +from spas import reconstruct +rec = reconstruct(model, device, meas_bin_2) ``` -Reconstruction and plot +Plot the spectral images ``` python -from spas import reconstruct, plot_color -rec = reconstruct(model, device, meas_bin[0:8192//4,:], 1, noise_bin) +from spas import plot_color plot_color(rec, wavelengths_bin) ``` diff --git a/scripts/reconstruction_neural_network.py b/scripts/reconstruction_neural_network.py index 9c1b786..5dd72ce 100644 --- a/scripts/reconstruction_neural_network.py +++ b/scripts/reconstruction_neural_network.py @@ -1,5 +1,11 @@ # -*- coding: utf-8 -*- """ +In this example, we consider the reconstruction of the 'Cat_LinearColoredFilter' acquisition that belongs to the SPIHIM collection. The acquisition was done at resolution 64x64 and the reconstrcution is performed at resolution 128x128. + +* The raw data can be downloaded [here](https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/622b5ea843258e76eab21740) + * The reconstruction model can be downloaded [here](https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/622b5ea843258e76eab21740) + * The reconstruction covariance matrix can be downloaded [here](?) + Created on Thu Jan 26 16:49:14 2023 @author: ducros @@ -50,7 +56,7 @@ _, acquisition_param, _, _ = read_metadata(meta_path) wavelengths = acquisition_param.wavelengths -# bin raw data +# bin raw data between 530 and 730 nm meas_bin, wavelengths_bin, _ = spectral_binning(meas.T, wavelengths, 530, 730, 4) #%% Reorder and subsample @@ -59,7 +65,7 @@ #%% Reconstruct from spas import reconstruct -rec = reconstruct(model, device, meas_bin_2, 1) +rec = reconstruct(model, device, meas_bin_2) #%% Plot from spas import plot_color diff --git a/spas/reconstruction_nn.py b/spas/reconstruction_nn.py index 81a3999..bc513df 100644 --- a/spas/reconstruction_nn.py +++ b/spas/reconstruction_nn.py @@ -29,24 +29,14 @@ from spyrit.core.Acquisition import Acquisition_Poisson_approx_Gauss from spyrit.core.Forward_Operator import Forward_operator_Split_ft_had from spyrit.core.Preprocess import Preprocess_Split_diag_poisson -from spyrit.core.Data_Consistency import Generalized_Orthogonal_Tikhonov, Pinv_orthogonal -from spyrit.core.neural_network import Unet, Identity +from spyrit.core.Data_Consistency import Generalized_Orthogonal_Tikhonov #, Pinv_orthogonal +from spyrit.core.neural_network import Unet #, Identity from spyrit.misc.sampling import Permutation_Matrix - from spas.noise import noiseClass from spas.metadata import AcquisitionParameters -class netType(IntEnum): - """Possible model architectures. - """ - c0mp = 0 - comp = 1 - pinv = 2 - free = 3 - - @dataclass class ReconstructionParameters: """Reconstruction parameters for loading a reconstruction model. @@ -148,9 +138,10 @@ def setup_reconstruction(cov_path: str, network_params.step_size, network_params.gamma, network_params.batch_size, network_params.regularization) - torch.cuda.empty_cache() # keep this here? + torch.cuda.empty_cache() # need to keep this here? title = Path(model_folder) / net_folder / (net_folder + suffix) load_net(title, model, device) + model.eval() # Mandantory when batchNorm is used model = model.to(device) return model, device @@ -158,6 +149,7 @@ def setup_reconstruction(cov_path: str, def reorder_subsample(meas: np.ndarray, acqui_param: AcquisitionParameters, recon_param: ReconstructionParameters, + recon_cov_path: str = "/path/cov.py", ) -> np.ndarray: """Reorder and subsample measurements @@ -166,9 +158,11 @@ def reorder_subsample(meas: np.ndarray, Spectral measurements with dimensions (N_wavelength x M_acq), where M_acq is the number of acquired patterns acqui_param (AcquisitionParameters): - Parameters used during the acquisition of the spectral measurments + Parameters used during the acquisition of the spectral measurements recon_param (ReconstructionParameters): Parameters of the reconstruction. + recon_cov_path (str, optional): + path to covariance matrix used for reconstruction Returns: (np.ndarray): @@ -194,12 +188,17 @@ def reorder_subsample(meas: np.ndarray, n_sub = math.ceil(recon_param.M**0.5) Ord_rec[:,n_sub:] = 0 Ord_rec[n_sub:,:] = 0 + + elif recon_param.subs == 'var': + Cov_rec = np.load(recon_cov_path) + Ord_rec = Cov2Var(Cov_rec) + Perm_rec = Permutation_Matrix(Ord_rec) # meas = meas.T - # Subsampled acquisition permutation matrix (and zero filling if necessary) + # Subsample acquisition permutation matrix (fill with zeros if necessary) if N_rec > N_acq: # Square subsampling in the "natural" order @@ -257,7 +256,7 @@ def reorder_subsample(meas: np.ndarray, def reconstruct(model: Union[Pinv_Net, DC2_Net], device: str, spectral_data: np.ndarray, - batches : int, + batches : int = 1, #noise_model: noiseClass = None, is_process: bool = False ) -> np.ndarray: @@ -274,26 +273,25 @@ def reconstruct(model: Union[Pinv_Net, DC2_Net], Device to which the model was loaded. Used to make sure the spectral data is loaded to the same device as the model. spectral_data (np.ndarray): - Spectral data acquired, must have the dimensions - (spectral dimension x patterns). - batches (int): + Spectral data acquired, must have the dimensions (spectral + dimension x patterns). + batches (int, optional): Number of batches for reconstruction in case the device does not - have enough memory to reconstruct all data at a single time. E.g. - when reconstructing data distributed along many wavelengths, there - may be too many data. - ?? noise_model (noiseClass, optional): - Loaded noise model in case the reconstruction should use a denoising - method. Defaults to None. + have enough memory to reconstruct all spectral channels at once. + Defaults to 1 (i.e., all spectral channels reconstructed at once). is_process (bool, optional): If True, reconstruction is performed in 'real-time' mode, using multiprocessing for reconstruction, thus allowing reconstruction and acquisition simultaneously. Defaults to False. Returns: - np.ndarray: Reconstructed images following the dimensions: - spectral dimension x image size x image size. + np.ndarray: Reconstructed images with dimensions (spectral dimension x + image size x image size). """ - + + # noise_model (noiseClass, optional): + # Loaded noise model in case the reconstruction should use a denoising + # method. Defaults to None. proportion = spectral_data.shape[0]//batches # Amount of wavelengths per batch @@ -337,9 +335,6 @@ def reconstruct(model: Union[Pinv_Net, DC2_Net], spectral_data_torch = torch.tensor(spectral_data[lambda_indices,:], dtype = torch.float, device = device) - - # spectral_data_torch = \ - # torch.Tensor(spectral_data[lambda_indices,:]).to(device) recon_torch = model.reconstruct_expe(spectral_data_torch)#, #len(lambda_indices), 1, model.n, model.n, C, s, K) From 419923e9f0d74c0add510e6094b8417a0046a949 Mon Sep 17 00:00:00 2001 From: nducros Date: Mon, 30 Jan 2023 18:41:43 +0100 Subject: [PATCH 03/11] minor update of doc --- README.md | 3 ++- scripts/reconstruction_neural_network.py | 18 +++++++++++++----- 2 files changed, 15 insertions(+), 6 deletions(-) diff --git a/README.md b/README.md index 11016f0..000fb16 100644 --- a/README.md +++ b/README.md @@ -269,7 +269,8 @@ model_folder = '../model/' model, device = setup_reconstruction(cov_path, model_folder, network_param) ``` -Load noise calibration parameters (provided with the data or computed using tools in `/noise-calibration`). :warning:: Not used anymore in the current implementation of `spas`. +Load noise calibration parameters (provided with the data or computed using tools in `/noise-calibration`). :warning: Noise parameters are not used anymore in the current implementation of `spas`. + ``` python from spas import load_noise noise = load_noise('../noise-calibration/fit_model.npz') diff --git a/scripts/reconstruction_neural_network.py b/scripts/reconstruction_neural_network.py index 5dd72ce..04b12af 100644 --- a/scripts/reconstruction_neural_network.py +++ b/scripts/reconstruction_neural_network.py @@ -1,11 +1,19 @@ # -*- coding: utf-8 -*- """ -In this example, we consider the reconstruction of the 'Cat_LinearColoredFilter' acquisition that belongs to the SPIHIM collection. The acquisition was done at resolution 64x64 and the reconstrcution is performed at resolution 128x128. +In this example, we consider the reconstruction of the 'Cat_LinearColoredFilter' +acquisition that belongs to the SPIHIM collection. The acquisition was done at +resolution 64x64 and the reconstrcution is performed at resolution 128x128. -* The raw data can be downloaded [here](https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/622b5ea843258e76eab21740) - * The reconstruction model can be downloaded [here](https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/622b5ea843258e76eab21740) - * The reconstruction covariance matrix can be downloaded [here](?) - +* The raw data can be downloaded here: + https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/622b5ea843258e76eab21740 + +* The reconstruction model can be downloaded here: + https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/638630794d15dd536f04831e + +* The covariance matrix can be downloaded here: + https://pilot-warehouse.creatis.insa-lyon.fr/#collection/6140ba6929e3fc10d47dbe3e/folder/63d7f3620386da2747641e1b + + Created on Thu Jan 26 16:49:14 2023 @author: ducros From ddf14b1918af6ea34788f9033d2ce2ff41650858 Mon Sep 17 00:00:00 2001 From: Mahieu-Williame Date: Fri, 31 Mar 2023 14:39:51 +0200 Subject: [PATCH 04/11] new figure name inside the overview folder ex: name_RGB_IMAGE_had_reco --- spas/metadata.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/spas/metadata.py b/spas/metadata.py index 0fa030a..ad38de1 100644 --- a/spas/metadata.py +++ b/spas/metadata.py @@ -1020,9 +1020,10 @@ def __init__(self, data_folder_name, data_name): self.data_path = self.subfolder_path + '/' + data_name self.had_reco_path = self.data_path + '_had_reco.npz' - self.fig_had_reco_path = self.overview_path + '/' + 'HAD_RECO_' + data_name + self.fig_had_reco_path = self.overview_path + '/' + data_name self.pathIDSsnapshot = Path(self.data_path + '_IDScam_before_acq.npy') self.pathIDSsnapshot_overview = self.overview_path + '/' + 'CAM_before_acq_' + data_name + '.png' - + self.nn_reco_path = self.data_path + '_nn_reco.npz' + self.fig_nn_reco_path = self.overview_path + '/' + data_name \ No newline at end of file From 0aba940b850fe32b6d0f777ac50b81b3d17649bd Mon Sep 17 00:00:00 2001 From: Mahieu-Williame Date: Fri, 31 Mar 2023 14:46:36 +0200 Subject: [PATCH 05/11] possibility or not to upload metadata --- spas/transfer_data_to_girder.py | 321 ++++++++++++++++---------------- 1 file changed, 164 insertions(+), 157 deletions(-) diff --git a/spas/transfer_data_to_girder.py b/spas/transfer_data_to_girder.py index 9eda1a3..456d54b 100644 --- a/spas/transfer_data_to_girder.py +++ b/spas/transfer_data_to_girder.py @@ -11,7 +11,7 @@ # from singlepixel import read_metadata def transfer_data(metadata, acquisition_parameters, spectrometer_params, DMD_params, - setup_version, data_folder_name, data_name): + setup_version, data_folder_name, data_name, upload_metadata): #%%########################## Girder info ################################# url = 'https://pilot-warehouse.creatis.insa-lyon.fr/api/v1' @@ -58,69 +58,75 @@ def transfer_data(metadata, acquisition_parameters, spectrometer_params, DMD_par #%%####################### prepare metadata dict ############################ #metadata, acquisition_parameters, spectrometer_params, DMD_params = read_metadata(data_path + '/' + data_name +'_metadata.json') - experiment_dict = metadata.__dict__ - acq_par_dict = acquisition_parameters.__dict__ - spectrometer_dict = spectrometer_params.__dict__ - DMD_params_dict = DMD_params.__dict__ - - experiment_dict2 = {} - for key in experiment_dict.keys(): - new_key = 'a)_EXP_' + key - experiment_dict2[new_key] = experiment_dict[key] + if upload_metadata == 1: + experiment_dict = metadata.__dict__ + acq_par_dict = acquisition_parameters.__dict__ + spectrometer_dict = spectrometer_params.__dict__ + DMD_params_dict = DMD_params.__dict__ - acq_par_dict2 = {} - for key in acq_par_dict.keys(): - new_key = 'b)_ACQ_' + key - #if not (key == 'patterns' or key == 'measurement_time' or key == 'timestamps' or key == 'wavelengths'): - acq_par_dict2[new_key] = acq_par_dict[key] - #print(key + '= ' + str(acq_par_dict[key])) - - spectrometer_dict2 = {} - for key in spectrometer_dict.keys(): - new_key = 'c)_SPECTRO_' + key - spectrometer_dict2[new_key] = spectrometer_dict[key] + experiment_dict2 = {} + for key in experiment_dict.keys(): + new_key = 'a)_EXP_' + key + experiment_dict2[new_key] = experiment_dict[key] + + acq_par_dict2 = {} + for key in acq_par_dict.keys(): + new_key = 'b)_ACQ_' + key + #if not (key == 'patterns' or key == 'measurement_time' or key == 'timestamps' or key == 'wavelengths'): + acq_par_dict2[new_key] = acq_par_dict[key] + #print(key + '= ' + str(acq_par_dict[key])) + + spectrometer_dict2 = {} + for key in spectrometer_dict.keys(): + new_key = 'c)_SPECTRO_' + key + spectrometer_dict2[new_key] = spectrometer_dict[key] + + DMD_params_dict2 = {} + for key in DMD_params_dict.keys(): + new_key = 'd)_DMD_' + key + DMD_params_dict2[new_key] = DMD_params_dict[key] + + dict = {} + dict.update(experiment_dict2) + dict.update(acq_par_dict2) + dict.update(spectrometer_dict2) + dict.update(DMD_params_dict2) - DMD_params_dict2 = {} - for key in DMD_params_dict.keys(): - new_key = 'd)_DMD_' + key - DMD_params_dict2[new_key] = DMD_params_dict[key] - - dict = {} - dict.update(experiment_dict2) - dict.update(acq_par_dict2) - dict.update(spectrometer_dict2) - dict.update(DMD_params_dict2) - - del dict['a)_EXP_output_directory'] - del dict['a)_EXP_pattern_order_source'] - del dict['a)_EXP_pattern_source'] - del dict['a)_EXP_class_description'] - del dict['b)_ACQ_class_description'] - del dict['b)_ACQ_patterns'] - del dict['b)_ACQ_measurement_time'] - del dict['b)_ACQ_timestamps'] - del dict['b)_ACQ_wavelengths'] - del dict['c)_SPECTRO_initial_available_pixels'] - del dict['c)_SPECTRO_store_to_ram'] - del dict['c)_SPECTRO_class_description'] - del dict['c)_SPECTRO_detector'] #SENS_HAMS11639 - del dict['c)_SPECTRO_firmware_version'] #001.011.000.000 - del dict['c)_SPECTRO_fpga_version'] #014.000.012.000 - del dict['c)_SPECTRO_dll_version'] #9.10.3.0 - del dict['d)_DMD_apps_fpga_temperature'] - del dict['d)_DMD_class_description'] - del dict['d)_DMD_ddc_fpga_temperature'] - del dict['d)_DMD_device_number'] - del dict['d)_DMD_id'] - del dict['d)_DMD_initial_memory'] - del dict['d)_DMD_pcb_temperature'] - del dict['d)_DMD_bitplanes'] - del dict['d)_DMD_type'] - del dict['d)_DMD_usb_connection'] - del dict['d)_DMD_ALP_version'] - - #%%################### begin metadata transfer ############################ - gc.addMetadataToFolder(data_folder_id, dict) + del dict['a)_EXP_output_directory'] + del dict['a)_EXP_pattern_order_source'] + del dict['a)_EXP_pattern_source'] + del dict['a)_EXP_class_description'] + del dict['b)_ACQ_class_description'] + del dict['b)_ACQ_patterns'] + try: + del dict['b)_ACQ_patterns_sub'] + except: + print('b)_ACQ_patterns_sub metadata does not exist') + del dict['b)_ACQ_patterns_wp'] + del dict['b)_ACQ_measurement_time'] + del dict['b)_ACQ_timestamps'] + del dict['b)_ACQ_wavelengths'] + del dict['c)_SPECTRO_initial_available_pixels'] + del dict['c)_SPECTRO_store_to_ram'] + del dict['c)_SPECTRO_class_description'] + del dict['c)_SPECTRO_detector'] #SENS_HAMS11639 + del dict['c)_SPECTRO_firmware_version'] #001.011.000.000 + del dict['c)_SPECTRO_fpga_version'] #014.000.012.000 + del dict['c)_SPECTRO_dll_version'] #9.10.3.0 + del dict['d)_DMD_apps_fpga_temperature'] + del dict['d)_DMD_class_description'] + del dict['d)_DMD_ddc_fpga_temperature'] + del dict['d)_DMD_device_number'] + del dict['d)_DMD_id'] + del dict['d)_DMD_initial_memory'] + del dict['d)_DMD_pcb_temperature'] + del dict['d)_DMD_bitplanes'] + del dict['d)_DMD_type'] + del dict['d)_DMD_usb_connection'] + del dict['d)_DMD_ALP_version'] + + #%%################### begin metadata transfer ############################ + gc.addMetadataToFolder(data_folder_id, dict) #%%##################### erase temp folder ################################ if len(os.listdir('../temp')) != 0: @@ -130,7 +136,7 @@ def transfer_data(metadata, acquisition_parameters, spectrometer_params, DMD_par def transfer_data_2arms(metadata, acquisition_parameters, spectrometer_params, DMD_params, camPar, - setup_version, data_folder_name, data_name): + setup_version, data_folder_name, data_name, upload_metadata): #unwrap structure into camPar camPar.AOI_X = camPar.rectAOI.s32X.value @@ -175,103 +181,104 @@ def transfer_data_2arms(metadata, acquisition_parameters, spectrometer_params, D #print('data_folder_id = ' + data_folder_id) #%%####################### prepare metadata dict ############################ #metadata, acquisition_parameters, spectrometer_params, DMD_params = read_metadata(data_path + '/' + data_name +'_metadata.json') - experiment_dict = metadata.__dict__ - acq_par_dict = acquisition_parameters.__dict__ - spectrometer_dict = spectrometer_params.__dict__ - DMD_params_dict = DMD_params.__dict__ - CAM_params_dict = camPar.__dict__ - - experiment_dict2 = {} - for key in experiment_dict.keys(): - new_key = 'a)_EXP_' + key - experiment_dict2[new_key] = experiment_dict[key] + if upload_metadata == 1: + experiment_dict = metadata.__dict__ + acq_par_dict = acquisition_parameters.__dict__ + spectrometer_dict = spectrometer_params.__dict__ + DMD_params_dict = DMD_params.__dict__ + CAM_params_dict = camPar.__dict__ - acq_par_dict2 = {} - for key in acq_par_dict.keys(): - new_key = 'b)_ACQ_' + key - #if not (key == 'patterns' or key == 'measurement_time' or key == 'timestamps' or key == 'wavelengths'): - acq_par_dict2[new_key] = acq_par_dict[key] - #print(key + '= ' + str(acq_par_dict[key])) - - spectrometer_dict2 = {} - for key in spectrometer_dict.keys(): - new_key = 'c)_SPECTRO_' + key - spectrometer_dict2[new_key] = spectrometer_dict[key] - - DMD_params_dict2 = {} - for key in DMD_params_dict.keys(): - new_key = 'd)_DMD_' + key - DMD_params_dict2[new_key] = DMD_params_dict[key] - - CAM_params_dict2 = {} - for key in CAM_params_dict.keys(): - if key == 'bandwidth': - new_key = 'e)_CAM_' + key + ' (MB/s)' - elif key == 'pixelClock': - new_key = 'e)_CAM_' + key + ' (MHz)' - elif key == 'exposureTime': - new_key = 'e)_CAM_' + key + ' (ms)' - else: - new_key = 'e)_CAM_' + key - CAM_params_dict2[new_key] = CAM_params_dict[key] + experiment_dict2 = {} + for key in experiment_dict.keys(): + new_key = 'a)_EXP_' + key + experiment_dict2[new_key] = experiment_dict[key] + + acq_par_dict2 = {} + for key in acq_par_dict.keys(): + new_key = 'b)_ACQ_' + key + #if not (key == 'patterns' or key == 'measurement_time' or key == 'timestamps' or key == 'wavelengths'): + acq_par_dict2[new_key] = acq_par_dict[key] + #print(key + '= ' + str(acq_par_dict[key])) + spectrometer_dict2 = {} + for key in spectrometer_dict.keys(): + new_key = 'c)_SPECTRO_' + key + spectrometer_dict2[new_key] = spectrometer_dict[key] - dict = {} - dict.update(experiment_dict2) - dict.update(acq_par_dict2) - dict.update(spectrometer_dict2) - dict.update(DMD_params_dict2) - dict.update(CAM_params_dict2) - - del dict['a)_EXP_output_directory'] - del dict['a)_EXP_pattern_order_source'] - del dict['a)_EXP_pattern_source'] - del dict['a)_EXP_class_description'] - del dict['b)_ACQ_class_description'] - del dict['b)_ACQ_patterns'] - del dict['b)_ACQ_patterns_wp'] - del dict['b)_ACQ_measurement_time'] - del dict['b)_ACQ_timestamps'] - del dict['b)_ACQ_wavelengths'] - del dict['c)_SPECTRO_initial_available_pixels'] - del dict['c)_SPECTRO_store_to_ram'] - del dict['c)_SPECTRO_class_description'] - del dict['c)_SPECTRO_detector']#SENS_HAMS11639 - del dict['c)_SPECTRO_firmware_version']#001.011.000.000 - del dict['c)_SPECTRO_fpga_version']#014.000.012.000 - del dict['c)_SPECTRO_dll_version']#9.10.3.0 - del dict['d)_DMD_apps_fpga_temperature'] - del dict['d)_DMD_class_description'] - del dict['d)_DMD_ddc_fpga_temperature'] - del dict['d)_DMD_device_number'] - del dict['d)_DMD_id'] - del dict['d)_DMD_initial_memory'] - del dict['d)_DMD_pcb_temperature'] - del dict['d)_DMD_bitplanes'] - del dict['d)_DMD_type'] - del dict['d)_DMD_usb_connection'] - del dict['d)_DMD_ALP_version'] - del dict['e)_CAM_hCam'] - del dict['e)_CAM_sInfo'] - del dict['e)_CAM_cInfo'] - del dict['e)_CAM_nBitsPerPixel'] - del dict['e)_CAM_m_nColorMode'] - del dict['e)_CAM_bytes_per_pixel'] - del dict['e)_CAM_rectAOI'] - del dict['e)_CAM_pcImageMemory'] - del dict['e)_CAM_MemID'] - del dict['e)_CAM_pitch'] - del dict['e)_CAM_camActivated'] - del dict['e)_CAM_Exit'] - del dict['e)_CAM_Memory'] - del dict['e)_CAM_avi'] - del dict['e)_CAM_punFileID'] - del dict['e)_CAM_timeout'] - del dict['e)_CAM_time_array'] - del dict['e)_CAM_black_pattern_num'] - - #%%################### begin metadata transfer ############################ - gc.addMetadataToFolder(data_folder_id, dict) + DMD_params_dict2 = {} + for key in DMD_params_dict.keys(): + new_key = 'd)_DMD_' + key + DMD_params_dict2[new_key] = DMD_params_dict[key] + + CAM_params_dict2 = {} + for key in CAM_params_dict.keys(): + if key == 'bandwidth': + new_key = 'e)_CAM_' + key + ' (MB/s)' + elif key == 'pixelClock': + new_key = 'e)_CAM_' + key + ' (MHz)' + elif key == 'exposureTime': + new_key = 'e)_CAM_' + key + ' (ms)' + else: + new_key = 'e)_CAM_' + key + CAM_params_dict2[new_key] = CAM_params_dict[key] + + + dict = {} + dict.update(experiment_dict2) + dict.update(acq_par_dict2) + dict.update(spectrometer_dict2) + dict.update(DMD_params_dict2) + dict.update(CAM_params_dict2) + + del dict['a)_EXP_output_directory'] + del dict['a)_EXP_pattern_order_source'] + del dict['a)_EXP_pattern_source'] + del dict['a)_EXP_class_description'] + del dict['b)_ACQ_class_description'] + del dict['b)_ACQ_patterns'] + del dict['b)_ACQ_patterns_wp'] + del dict['b)_ACQ_measurement_time'] + del dict['b)_ACQ_timestamps'] + del dict['b)_ACQ_wavelengths'] + del dict['c)_SPECTRO_initial_available_pixels'] + del dict['c)_SPECTRO_store_to_ram'] + del dict['c)_SPECTRO_class_description'] + del dict['c)_SPECTRO_detector']#SENS_HAMS11639 + del dict['c)_SPECTRO_firmware_version']#001.011.000.000 + del dict['c)_SPECTRO_fpga_version']#014.000.012.000 + del dict['c)_SPECTRO_dll_version']#9.10.3.0 + del dict['d)_DMD_apps_fpga_temperature'] + del dict['d)_DMD_class_description'] + del dict['d)_DMD_ddc_fpga_temperature'] + del dict['d)_DMD_device_number'] + del dict['d)_DMD_id'] + del dict['d)_DMD_initial_memory'] + del dict['d)_DMD_pcb_temperature'] + del dict['d)_DMD_bitplanes'] + del dict['d)_DMD_type'] + del dict['d)_DMD_usb_connection'] + del dict['d)_DMD_ALP_version'] + del dict['e)_CAM_hCam'] + del dict['e)_CAM_sInfo'] + del dict['e)_CAM_cInfo'] + del dict['e)_CAM_nBitsPerPixel'] + del dict['e)_CAM_m_nColorMode'] + del dict['e)_CAM_bytes_per_pixel'] + del dict['e)_CAM_rectAOI'] + del dict['e)_CAM_pcImageMemory'] + del dict['e)_CAM_MemID'] + del dict['e)_CAM_pitch'] + del dict['e)_CAM_camActivated'] + del dict['e)_CAM_Exit'] + del dict['e)_CAM_Memory'] + del dict['e)_CAM_avi'] + del dict['e)_CAM_punFileID'] + del dict['e)_CAM_timeout'] + del dict['e)_CAM_time_array'] + del dict['e)_CAM_black_pattern_num'] + + #%%################### begin metadata transfer ############################ + gc.addMetadataToFolder(data_folder_id, dict) #%%##################### erase temp folder ################################ if len(os.listdir('../temp')) != 0: list_TempFolder = os.listdir('../temp') From da03f20b777f559e80ea6c846c34cb91a5c2eba8 Mon Sep 17 00:00:00 2001 From: Mahieu-Williame Date: Fri, 31 Mar 2023 14:57:03 +0200 Subject: [PATCH 06/11] change plot_without_NN and with_NN functions --- spas/visualization.py | 114 ++++++++++++++++++++++++++++-------------- 1 file changed, 76 insertions(+), 38 deletions(-) diff --git a/spas/visualization.py b/spas/visualization.py index 6e05de7..adc0456 100644 --- a/spas/visualization.py +++ b/spas/visualization.py @@ -12,7 +12,7 @@ from matplotlib.colors import ListedColormap from matplotlib import pyplot as plt from mpl_toolkits.axes_grid1 import make_axes_locatable -from spas import plot_spec_to_rgb_image as plt_rgb +from spas.plot_spec_to_rgb_image import plot_spec_to_rgb_image from .noise import noiseClass # Libraries for the IDS CAMERA @@ -183,8 +183,8 @@ def spectral_slicing(F: np.ndarray, wavelengths: np.ndarray, lambda_min: int, of mu, sigma and K according to the generated wavelength bins. """ - if F.ndim == 3: - F = np.transpose(F,(0, 2, 1)) + # if F.ndim == 3: + # F = np.transpose(F,(0, 2, 1)) # Find the indices that fit the spectral range lambda_ind, = np.where((wavelengths > lambda_min) & @@ -475,49 +475,52 @@ def displayVid(camPar): break -def plot_reco_without_NN(acquisition_parameters, GT, Q, had_reco_path, fig_had_reco_path): +def plot_reco_without_NN(acquisition_parameters, GT, Q, all_path): + + had_reco_path = all_path.had_reco_path + fig_had_reco_path = all_path.fig_had_reco_path - GT = np.rot90(GT, 1) - GT = np.rot90(GT, 1) - if not os.path.exists(had_reco_path): np.savez_compressed(had_reco_path, GT) + + GT = np.rot90(GT, 1) + GT = np.rot90(GT, 1) F_bin, wavelengths_bin, bin_width = spectral_binning(GT.T, acquisition_parameters.wavelengths, 530, 730, 8) F_bin_rot = np.rot90(F_bin, axes=(1,2)) F_bin_flip = F_bin_rot[:,::-1,:] F_bin_1px, wavelengths_bin, bin_width = spectral_slicing(GT.T, acquisition_parameters.wavelengths, 530, 730, 8) - + F_bin_1px_rot = np.rot90(F_bin_1px, axes=(1,2)) + F_bin_1px_flip = F_bin_1px_rot[:,::-1,:] ############### spatial view, wavelength bin ############# #plt.figure() plot_color(F_bin_flip, wavelengths_bin) - plt.savefig(fig_had_reco_path + '_spatial_view_sum_wavelength_binning.png') + plt.savefig(fig_had_reco_path + '_BIN_IMAGE_had_reco.png') plt.show() ############### spatial view, one wavelength ############# #plt.figure() - plot_color(F_bin_1px, wavelengths_bin) - plt.savefig(fig_had_reco_path + '_spatial_view_single_slide_by_wavelength.png') + plot_color(F_bin_1px_flip, wavelengths_bin) + plt.savefig(fig_had_reco_path + '_SLICE_IMAGE_had_reco.png') plt.show() ############### spatial view, wavelength sum ############# #plt.figure() plt.imshow(np.sum(GT[:,:,193:877], axis=2))#[:,:,193:877] #(540-625 nm) plt.title('Sum of all wavelengths') - plt.savefig(fig_had_reco_path + '_spatial_view_sum_of_wavelengths.png') + plt.savefig(fig_had_reco_path + '_GRAY_IMAGE_had_reco.png') plt.show() ####################### RGB view ######################## print('Beging RGB convertion ...') - image_arr = plt_rgb.plot_spec_to_rgb_image(GT, acquisition_parameters.wavelengths) + image_arr = plot_spec_to_rgb_image(GT, acquisition_parameters.wavelengths) print('RGB convertion finished') plt.figure() - plt.imshow(image_arr, extent=[0, 10.5, 0, 10.5]) - plt.xlabel('X (mm)') - plt.ylabel('Y (mm)') - plt.savefig(fig_had_reco_path + '_RGB_view.png') + plt.imshow(image_arr) #, extent=[0, 10.5, 0, 10.5]) + # plt.xlabel('X (mm)') + # plt.ylabel('Y (mm)') + plt.savefig(fig_had_reco_path + '_RGB_IMAGE_had_reco.png') plt.show() - ####################### spectral view ################### GT50 = GT[16:48,16:48,:] GT25 = GT[24:40,24:40,:] @@ -529,31 +532,66 @@ def plot_reco_without_NN(acquisition_parameters, GT, Q, had_reco_path, fig_had_r plt.title("% of region from the center of the image") plt.legend(['25%', '50%', '100%']) plt.xlabel(r'$\lambda$ (nm)') - plt.savefig(fig_had_reco_path + '_spectral_view.png') + plt.savefig(fig_had_reco_path + '_SPECTRA_PLOT_had_reco.png') plt.show() -def plot_reco_with_NN(acquisition_parameters, network_params, spectral_data, noise, model, device): - CRreco = int(acquisition_parameters.acquired_spectra / 2 / network_params.CR) - F_bin, wavelengths_bin, bin_width, noise_bin = spectral_binning(spectral_data.T, acquisition_parameters.wavelengths, 530, 730, 8, noise) - recon = reconstruct(model, device, F_bin[:,0:acquisition_parameters.acquired_spectra//CRreco], 4, noise_bin) - # recon = reconstruct(model, device, F_bin[:,0:8192//4], 4, noise_bin) - recon = np.flip(recon, axis = 1) - recon = np.flip(recon, axis = 2) - plot_color(recon, wavelengths_bin) +def plot_reco_with_NN(acquisition_parameters, reco, all_path): + + reco = reco.transpose(1, 2, 0) + + nn_reco_path = all_path.nn_reco_path + fig_nn_reco_path = all_path.fig_nn_reco_path + + if not os.path.exists(nn_reco_path): + np.savez_compressed(nn_reco_path, reco) + + reco = np.rot90(reco, 1) + reco = np.rot90(reco, 1) + + F_bin, wavelengths_bin, bin_width = spectral_binning(reco.T, acquisition_parameters.wavelengths, 530, 730, 8) + F_bin_rot = np.rot90(F_bin, axes=(1,2)) + F_bin_flip = F_bin_rot[:,::-1,:] + F_bin_1px, wavelengths_bin, bin_width = spectral_slicing(reco.T, acquisition_parameters.wavelengths, 530, 730, 8) + F_bin_1px_rot = np.rot90(F_bin_1px, axes=(1,2)) + F_bin_1px_flip = F_bin_1px_rot[:,::-1,:] + ############### spatial view, wavelength bin ############# + #plt.figure() + plot_color(F_bin_flip, wavelengths_bin) + plt.savefig(fig_nn_reco_path + '_BIN_IMAGE_nn_reco.png') + plt.show() + + ############### spatial view, one wavelength ############# + #plt.figure() + plot_color(F_bin_1px_flip, wavelengths_bin) + plt.savefig(fig_nn_reco_path + '_SLICE_IMAGE_nn_reco.png') plt.show() - plt.imshow(np.sum(recon, axis=0)) - plt.title('NN reco, sum of all wavelengths') + ############### spatial view, wavelength sum ############# + #plt.figure() + plt.imshow(np.sum(reco[:,:,193:877], axis=2))#[:,:,193:877] #(540-625 nm) + plt.title('Sum of all wavelengths') + plt.savefig(fig_nn_reco_path + '_GRAY_IMAGE_nn_reco.png') plt.show() - F_bin, wavelengths_bin, bin_width, noise_bin = spectral_slicing(spectral_data.T, acquisition_parameters.wavelengths, 530, 730, 8, noise) - recon2 = reconstruct(model, device, F_bin[:,0:acquisition_parameters.acquired_spectra//CRreco], 4, noise_bin) - # recon2 = reconstruct(model, device, F_bin[:,0:8192//4], 4, noise_bin) - recon2 = np.flip(recon2, axis = 1) - recon2 = np.flip(recon2, axis = 2) - plot_color(recon2, wavelengths_bin) + ####################### RGB view ######################## + print('Beging RGB convertion ...') + image_arr = plot_spec_to_rgb_image(reco, acquisition_parameters.wavelengths) + print('RGB convertion finished') + plt.figure() + plt.imshow(image_arr) + plt.savefig(fig_nn_reco_path + '_RGB_IMAGE_nn_reco.png') plt.show() - - # if save_path: - # fig.savefig(f'{save_path}.png', dpi=300, bbox_inches='tight') + ####################### spectral view ################### + GT50 = reco[16:48,16:48,:] + GT25 = reco[24:40,24:40,:] + plt.figure() + plt.plot(acquisition_parameters.wavelengths, np.mean(np.mean(GT25,axis=1),axis=0)) + plt.plot(acquisition_parameters.wavelengths, np.mean(np.mean(GT50,axis=1),axis=0)) + plt.plot(acquisition_parameters.wavelengths, np.mean(np.mean(reco,axis=1),axis=0)) + plt.grid() + plt.title("% of region from the center of the image") + plt.legend(['25%', '50%', '100%']) + plt.xlabel(r'$\lambda$ (nm)') + plt.savefig(fig_nn_reco_path + '_SPECTRA_PLOT_nn_reco.png') + plt.show() \ No newline at end of file From 0c6391905644c48e1c63148bfef8caf96cc73808 Mon Sep 17 00:00:00 2001 From: Mahieu-Williame Date: Mon, 24 Apr 2023 13:51:33 +0200 Subject: [PATCH 07/11] Make it compatible with spyritv2 --- spas/reconstruction_nn.py | 150 ++++++++++++++++++++++++++------------ 1 file changed, 102 insertions(+), 48 deletions(-) diff --git a/spas/reconstruction_nn.py b/spas/reconstruction_nn.py index bc513df..b5d3242 100644 --- a/spas/reconstruction_nn.py +++ b/spas/reconstruction_nn.py @@ -20,18 +20,32 @@ import numpy as np from matplotlib import pyplot as plt -from spyrit.core.reconstruction import Pinv_Net, DC2_Net +# from spyrit.core.reconstruction import Pinv_Net, DC2_Net +from spyrit.core.recon import PinvNet, DCNet # modified by LMW 30/03/2023 + +# from spyrit.core.training import load_net +from spyrit.core.train import load_net # modified by LMW 30/03/2023 -from spyrit.core.training import load_net from spyrit.misc.statistics import Cov2Var from spyrit.misc.walsh_hadamard import walsh2_matrix -from spyrit.core.Acquisition import Acquisition_Poisson_approx_Gauss -from spyrit.core.Forward_Operator import Forward_operator_Split_ft_had -from spyrit.core.Preprocess import Preprocess_Split_diag_poisson -from spyrit.core.Data_Consistency import Generalized_Orthogonal_Tikhonov #, Pinv_orthogonal -from spyrit.core.neural_network import Unet #, Identity -from spyrit.misc.sampling import Permutation_Matrix +# from spyrit.core.Acquisition import Acquisition_Poisson_approx_Gauss +from spyrit.core.noise import Poisson # modified by LMW 30/03/2023 + +# from spyrit.core.Forward_Operator import Forward_operator_Split_ft_had +from spyrit.core.meas import HadamSplit # modified by LMW 30/03/2023 + +# from spyrit.core.Preprocess import Preprocess_Split_diag_poisson +from spyrit.core.prep import SplitPoisson # modified by LMW 30/03/2023 +# il y a aussi la classe SplitRowPoisson ??? + +# from spyrit.core.Data_Consistency import Generalized_Orthogonal_Tikhonov #, Pinv_orthogonal +from spyrit.core.recon import TikhonovMeasurementPriorDiag # modified by LMW 30/03/2023 + +# from spyrit.core.neural_network import Unet #, Identity +from spyrit.core.nnet import Unet # modified by LMW 30/03/2023 + +from spyrit.misc.sampling import Permutation_Matrix, reorder from spas.noise import noiseClass @@ -67,7 +81,7 @@ class ReconstructionParameters: def setup_reconstruction(cov_path: str, model_folder: str, network_params: ReconstructionParameters - ) -> Tuple[Union[Pinv_Net,DC2_Net], str]: + ) -> Tuple[Union[PinvNet,DCNet], str]: """Loads a neural network for reconstruction. Limited to measurements from patterns of size 2**K for reconstruction at @@ -106,53 +120,90 @@ def setup_reconstruction(cov_path: str, Ord[n_sub:,:] = 0 # Init network - Perm_rec = Permutation_Matrix(Ord) - Hperm = Perm_rec @ H - Pmat = Hperm[:network_params.M,:] + #Perm_rec = Permutation_Matrix(Ord) + #Hperm = Perm_rec @ H + #Pmat = Hperm[:network_params.M,:] # init - Forward = Forward_operator_Split_ft_had(Pmat, Perm_rec, - network_params.img_size, - network_params.img_size) - Noise = Acquisition_Poisson_approx_Gauss(network_params.N0, Forward) - Prep = Preprocess_Split_diag_poisson(network_params.N0, + # Forward = Forward_operator_Split_ft_had(Pmat, Perm_rec, + # network_params.img_size, + # network_params.img_size) + + Forward = HadamSplit(network_params.M, + network_params.img_size, + Ord)# modified by LMW 30/03/2023 + + # Noise = Acquisition_Poisson_approx_Gauss(network_params.N0, Forward) + + Noise = Poisson(Forward, network_params.N0) + + # Prep = Preprocess_Split_diag_poisson(network_params.N0, + # network_params.M, + # network_params.img_size**2) + + Prep = SplitPoisson(network_params.N0, network_params.M, network_params.img_size**2) Denoi = Unet() - Cov_perm = Perm_rec @ Cov_rec @ Perm_rec.T - DC = Generalized_Orthogonal_Tikhonov(sigma_prior = Cov_perm, - M = network_params.M, - N = network_params.img_size**2) - model = DC2_Net(Noise, Prep, DC, Denoi) - - # load - net_folder = '{}_{}_{}'.format( - network_params.arch, network_params.denoi, - network_params.data) - - suffix = '_{}_N0_{}_N_{}_M_{}_epo_{}_lr_{}_sss_{}_sdr_{}_bs_{}_reg_{}'.format( - network_params.subs, network_params.N0, - network_params.img_size, network_params.M, - network_params.num_epochs, network_params.learning_rate, - network_params.step_size, network_params.gamma, - network_params.batch_size, network_params.regularization) - - torch.cuda.empty_cache() # need to keep this here? - title = Path(model_folder) / net_folder / (net_folder + suffix) - load_net(title, model, device) - model.eval() # Mandantory when batchNorm is used - model = model.to(device) + # Cov_perm = Perm_rec @ Cov_rec @ Perm_rec.T + # DC = Generalized_Orthogonal_Tikhonov(sigma_prior = Cov_perm, + # M = network_params.M, + # N = network_params.img_size**2) + + + + # model = DC2_Net(Noise, Prep, DC, Denoi) + model = DCNet(Noise, Prep, Cov_rec, Denoi) + + # # load + # net_folder = '{}_{}_{}'.format( + # network_params.arch, network_params.denoi, + # network_params.data) + + # suffix = '_{}_N0_{}_N_{}_M_{}_epo_{}_lr_{}_sss_{}_sdr_{}_bs_{}_reg_{}'.format( + # network_params.subs, network_params.N0, + # network_params.img_size, network_params.M, + # network_params.num_epochs, network_params.learning_rate, + # network_params.step_size, network_params.gamma, + # network_params.batch_size, network_params.regularization) + + # torch.cuda.empty_cache() # need to keep this here? + # title = Path(model_folder) / net_folder / (net_folder + suffix) + # load_net(title, model, device) + # model.eval() # Mandantory when batchNorm is used + # model = model.to(device) + + # Load trained DC-Net + net_arch = network_params.arch + net_denoi = network_params.denoi + net_data = network_params.data + if (network_params.img_size == 128) and (network_params.M == 4096): + net_order = 'rect' + else: + net_order = 'var' + + bs = 256 + # net_suffix = f'N0_{network_params.N0}_N_{network_params.img_size}_M_{network_params.M}_epo_30_lr_0.001_sss_10_sdr_0.5_bs_{bs}_reg_1e-07_light' + net_suffix = f'N0_{network_params.N0}_N_{network_params.img_size}_M_{network_params.M}_epo_30_lr_0.001_sss_10_sdr_0.5_bs_{bs}_reg_1e-07_light' + + net_folder= f'{net_arch}_{net_denoi}_{net_data}/' + + net_title = f'{net_arch}_{net_denoi}_{net_data}_{net_order}_{net_suffix}' + # title = './model_v2/' + net_folder + net_title + title = 'C:/openspyrit/models/' + net_folder + net_title + load_net(title, model, device, False) + model.eval() # Mandantory when batchNorm is used return model, device + def reorder_subsample(meas: np.ndarray, acqui_param: AcquisitionParameters, recon_param: ReconstructionParameters, recon_cov_path: str = "/path/cov.py", ) -> np.ndarray: """Reorder and subsample measurements - Args: meas (np.ndarray): Spectral measurements with dimensions (N_wavelength x M_acq), where @@ -163,13 +214,11 @@ def reorder_subsample(meas: np.ndarray, Parameters of the reconstruction. recon_cov_path (str, optional): path to covariance matrix used for reconstruction - Returns: (np.ndarray): Spectral measurements with dimensions (N_wavelength x M_rec), where M_rec is the number of patterns considered for reconstruction. Acquisitions can be subsampled a posteriori, leadind to M_rec < M_acq - """ # Dimensions (N.B: images are assumed to be square) N_acq = acqui_param.pattern_dimension_x @@ -253,7 +302,7 @@ def reorder_subsample(meas: np.ndarray, return meas[:2*recon_param.M,:].T -def reconstruct(model: Union[Pinv_Net, DC2_Net], +def reconstruct(model: Union[PinvNet, DCNet], device: str, spectral_data: np.ndarray, batches : int = 1, @@ -295,13 +344,18 @@ def reconstruct(model: Union[Pinv_Net, DC2_Net], proportion = spectral_data.shape[0]//batches # Amount of wavelengths per batch - img_size = model.Acq.FO.h # image assumed to be square + # img_size = model.Acq.FO.h # image assumed to be square + # img_size = 128 # Modified by LMW 30/03/2023 + img_size = model.Acq.meas_op.h + recon = np.zeros((spectral_data.shape[0], img_size, img_size)) start = perf_counter_ns() - model.PreP.set_expe() - + # model.PreP.set_expe() + model.prep.set_expe() # Modified by LMW 30/03/2023 + model.to(device) # Modified by LMW 30/03/2023 + with torch.no_grad(): for batch in range(batches): @@ -353,7 +407,7 @@ def reconstruct(model: Union[Pinv_Net, DC2_Net], return recon -def reconstruct_process(model: Union[Pinv_Net, DC2_Net], +def reconstruct_process(model: Union[PinvNet, DCNet], device: str, queue_to_recon: Queue, queue_reconstructed: Queue, batches: int, noise_model: noiseClass, sleep_time: float = 0.3) -> None: """Performs reconstruction in real-time. From d38a8a896ce0b272d2e51797e334d93516fd18fc Mon Sep 17 00:00:00 2001 From: nducros Date: Thu, 27 Apr 2023 17:20:42 +0200 Subject: [PATCH 08/11] compatible with spyrit 2.1 --- spas/reconstruction_nn.py | 150 ++++---------------------------------- 1 file changed, 15 insertions(+), 135 deletions(-) diff --git a/spas/reconstruction_nn.py b/spas/reconstruction_nn.py index b5d3242..37035a5 100644 --- a/spas/reconstruction_nn.py +++ b/spas/reconstruction_nn.py @@ -20,34 +20,16 @@ import numpy as np from matplotlib import pyplot as plt -# from spyrit.core.reconstruction import Pinv_Net, DC2_Net -from spyrit.core.recon import PinvNet, DCNet # modified by LMW 30/03/2023 - -# from spyrit.core.training import load_net -from spyrit.core.train import load_net # modified by LMW 30/03/2023 - +from spyrit.core.train import load_net +from spyrit.core.noise import Poisson +from spyrit.core.meas import HadamSplit +from spyrit.core.prep import SplitPoisson +from spyrit.core.recon import PinvNet, DCNet, TikhonovMeasurementPriorDiag +from spyrit.core.nnet import Unet +from spyrit.misc.sampling import Permutation_Matrix, reorder from spyrit.misc.statistics import Cov2Var from spyrit.misc.walsh_hadamard import walsh2_matrix -# from spyrit.core.Acquisition import Acquisition_Poisson_approx_Gauss -from spyrit.core.noise import Poisson # modified by LMW 30/03/2023 - -# from spyrit.core.Forward_Operator import Forward_operator_Split_ft_had -from spyrit.core.meas import HadamSplit # modified by LMW 30/03/2023 - -# from spyrit.core.Preprocess import Preprocess_Split_diag_poisson -from spyrit.core.prep import SplitPoisson # modified by LMW 30/03/2023 -# il y a aussi la classe SplitRowPoisson ??? - -# from spyrit.core.Data_Consistency import Generalized_Orthogonal_Tikhonov #, Pinv_orthogonal -from spyrit.core.recon import TikhonovMeasurementPriorDiag # modified by LMW 30/03/2023 - -# from spyrit.core.neural_network import Unet #, Identity -from spyrit.core.nnet import Unet # modified by LMW 30/03/2023 - -from spyrit.misc.sampling import Permutation_Matrix, reorder - - from spas.noise import noiseClass from spas.metadata import AcquisitionParameters @@ -119,60 +101,17 @@ def setup_reconstruction(cov_path: str, Ord[:,n_sub:] = 0 Ord[n_sub:,:] = 0 - # Init network - #Perm_rec = Permutation_Matrix(Ord) - #Hperm = Perm_rec @ H - #Pmat = Hperm[:network_params.M,:] - - # init - # Forward = Forward_operator_Split_ft_had(Pmat, Perm_rec, - # network_params.img_size, - # network_params.img_size) - + # Init network Forward = HadamSplit(network_params.M, network_params.img_size, - Ord)# modified by LMW 30/03/2023 - - # Noise = Acquisition_Poisson_approx_Gauss(network_params.N0, Forward) - + Ord) Noise = Poisson(Forward, network_params.N0) - - # Prep = Preprocess_Split_diag_poisson(network_params.N0, - # network_params.M, - # network_params.img_size**2) - Prep = SplitPoisson(network_params.N0, - network_params.M, - network_params.img_size**2) - + network_params.M, + network_params.img_size**2) Denoi = Unet() - # Cov_perm = Perm_rec @ Cov_rec @ Perm_rec.T - # DC = Generalized_Orthogonal_Tikhonov(sigma_prior = Cov_perm, - # M = network_params.M, - # N = network_params.img_size**2) - - - - # model = DC2_Net(Noise, Prep, DC, Denoi) model = DCNet(Noise, Prep, Cov_rec, Denoi) - # # load - # net_folder = '{}_{}_{}'.format( - # network_params.arch, network_params.denoi, - # network_params.data) - - # suffix = '_{}_N0_{}_N_{}_M_{}_epo_{}_lr_{}_sss_{}_sdr_{}_bs_{}_reg_{}'.format( - # network_params.subs, network_params.N0, - # network_params.img_size, network_params.M, - # network_params.num_epochs, network_params.learning_rate, - # network_params.step_size, network_params.gamma, - # network_params.batch_size, network_params.regularization) - - # torch.cuda.empty_cache() # need to keep this here? - # title = Path(model_folder) / net_folder / (net_folder + suffix) - # load_net(title, model, device) - # model.eval() # Mandantory when batchNorm is used - # model = model.to(device) # Load trained DC-Net net_arch = network_params.arch @@ -184,13 +123,11 @@ def setup_reconstruction(cov_path: str, net_order = 'var' bs = 256 - # net_suffix = f'N0_{network_params.N0}_N_{network_params.img_size}_M_{network_params.M}_epo_30_lr_0.001_sss_10_sdr_0.5_bs_{bs}_reg_1e-07_light' net_suffix = f'N0_{network_params.N0}_N_{network_params.img_size}_M_{network_params.M}_epo_30_lr_0.001_sss_10_sdr_0.5_bs_{bs}_reg_1e-07_light' net_folder= f'{net_arch}_{net_denoi}_{net_data}/' net_title = f'{net_arch}_{net_denoi}_{net_data}_{net_order}_{net_suffix}' - # title = './model_v2/' + net_folder + net_title title = 'C:/openspyrit/models/' + net_folder + net_title load_net(title, model, device, False) model.eval() # Mandantory when batchNorm is used @@ -230,7 +167,6 @@ def reorder_subsample(meas: np.ndarray, Ord_acq = np.reshape(Ord_acq, (N_acq,N_acq)) # sampling map Perm_acq = Permutation_Matrix(Ord_acq).T - # Order used for reconstruction if recon_param.subs == 'rect': Ord_rec = np.ones((N_rec, N_rec)) @@ -244,61 +180,10 @@ def reorder_subsample(meas: np.ndarray, Perm_rec = Permutation_Matrix(Ord_rec) - # + # reorder meas = meas.T - - # Subsample acquisition permutation matrix (fill with zeros if necessary) - if N_rec > N_acq: - - # Square subsampling in the "natural" order - Ord_sub = np.zeros((N_rec,N_rec)) - Ord_sub[:N_acq,:N_acq]= -np.arange(-N_acq**2,0).reshape(N_acq,N_acq) - Perm_sub = Permutation_Matrix(Ord_sub) - - # Natural order measurements (N_acq resolution) - Perm_raw = np.zeros((2*N_acq**2,2*N_acq**2)) - Perm_raw[::2,::2] = Perm_acq.T - Perm_raw[1::2,1::2] = Perm_acq.T - meas = Perm_raw @ meas - - # Zero filling (needed only when reconstruction resolution is higher - # than acquisition res) - zero_filled = np.zeros((2*N_rec**2, N_wav)) - zero_filled[:2*N_acq**2,:] = meas - - meas = zero_filled + meas = reorder(meas, Perm_acq, Perm_rec) - Perm_raw = np.zeros((2*N_rec**2,2*N_rec**2)) - Perm_raw[::2,::2] = Perm_sub.T - Perm_raw[1::2,1::2] = Perm_sub.T - - meas = Perm_raw @ meas - - elif N_rec == N_acq: - Perm_sub = Perm_acq[:N_rec**2,:].T - - elif N_rec < N_acq: - # Square subsampling in the "natural" order - Ord_sub = np.zeros((N_acq,N_acq)) - Ord_sub[:N_rec,:N_rec]= -np.arange(-N_rec**2,0).reshape(N_rec,N_rec) - Perm_sub = Permutation_Matrix(Ord_sub) - Perm_sub = Perm_sub[:N_rec**2,:] - Perm_sub = Perm_sub @ Perm_acq.T - - #Reorder measurements when reconstruction order is not "natural" - if N_rec <= N_acq: - # Get both positive and negative coefficients permutated - Perm = Perm_rec @ Perm_sub - Perm_raw = np.zeros((2*N_rec**2,2*N_acq**2)) - - elif N_rec > N_acq: - Perm = Perm_rec - Perm_raw = np.zeros((2*N_rec**2,2*N_rec**2)) - - Perm_raw[::2,::2] = Perm - Perm_raw[1::2,1::2] = Perm - meas = Perm_raw @ meas - return meas[:2*recon_param.M,:].T @@ -343,18 +228,13 @@ def reconstruct(model: Union[PinvNet, DCNet], # method. Defaults to None. proportion = spectral_data.shape[0]//batches # Amount of wavelengths per batch - - # img_size = model.Acq.FO.h # image assumed to be square - # img_size = 128 # Modified by LMW 30/03/2023 img_size = model.Acq.meas_op.h - recon = np.zeros((spectral_data.shape[0], img_size, img_size)) - start = perf_counter_ns() # model.PreP.set_expe() - model.prep.set_expe() # Modified by LMW 30/03/2023 - model.to(device) # Modified by LMW 30/03/2023 + model.prep.set_expe() + model.to(device) with torch.no_grad(): for batch in range(batches): From 2672d0e51c2bf5c4410d72591417a1cb1fc667ec Mon Sep 17 00:00:00 2001 From: Mahieu-Williame Date: Wed, 10 May 2023 12:17:10 +0200 Subject: [PATCH 09/11] DLLs are placed in the function : TRY, to be used by Linux and/or MacOS in the case the DLL are not installed --- scripts/main_seq_2arms.py | 127 ++-- ...reconstruct_and_transfer_data_to_girder.py | 577 ++++++++++++++++++ spas/acquisition.py | 31 +- spas/metadata.py | 21 +- spas/reconstruction_nn.py | 37 +- spas/transfer_data_to_girder.py | 14 +- spas/visualization.py | 81 ++- 7 files changed, 750 insertions(+), 138 deletions(-) create mode 100644 scripts/reconstruct_and_transfer_data_to_girder.py diff --git a/scripts/main_seq_2arms.py b/scripts/main_seq_2arms.py index 4cac646..405736e 100644 --- a/scripts/main_seq_2arms.py +++ b/scripts/main_seq_2arms.py @@ -9,19 +9,21 @@ from spas.acquisition import init_2arms, setup_cam, AcquisitionParameters, setup_2arms, acquire, acquire_2arms, snapshot, disconnect_2arms, captureVid from spas.metadata import MetaData, func_path from spas.reconstruction import reconstruction_hadamard -from spas.reconstruction_nn import ReconstructionParameters, setup_reconstruction +from spas.reconstruction_nn import ReconstructionParameters, setup_reconstruction, reorder_subsample from spas.noise import load_noise from spas.visualization import snapshotVisu, displayVid, plot_reco_without_NN, plot_reco_with_NN from spas.transfer_data_to_girder import transfer_data_2arms import spyrit.misc.walsh_hadamard as wh +from spas import reconstruct +import time #%% Initialize hardware spectrometer, DMD, DMD_initial_memory, camPar = init_2arms() #%% Define the AOI # Warning, not all values are allowed for Width and Height (max: 2076x3088 | ex: 768x544) -camPar.rectAOI.s32X.value = 880# // X -camPar.rectAOI.s32Y.value = 590# // Y -camPar.rectAOI.s32Width.value = 768#1544 # 3088 // Width must be multiple of 8 -camPar.rectAOI.s32Height.value = 544#1038# 2076# // Height +camPar.rectAOI.s32X.value = 0 #0#1100# // X +camPar.rectAOI.s32Y.value = 0#0#800#765# // Y +camPar.rectAOI.s32Width.value = 1544 #3000#3088#3088# 768 // Width must be multiple of 8 +camPar.rectAOI.s32Height.value = 1730#2000## 1038# 544 # 2076 // Height camPar = captureVid(camPar) #%% Set Camera Parameters @@ -32,41 +34,44 @@ Gain = 0, # Gain boundary : [0 100] gain_boost = 'OFF', # set "ON" to activate gain boost, "OFF" to deactivate nGamma = 1, # Gamma boundary : [1 - 2.2] - ExposureTime = 0.09, # Exposure Time (ms) boudary : [0.013 - 56.221] - black_level = 0) # Black Level boundary : [0 255] + ExposureTime = 1,# Exposure Time (ms) boudary : [0.013 - 56.221] + black_level = 5) # Black Level boundary : [0 255] snapshotVisu(camPar) #%% Display video in continous mode for optical tuning displayVid(camPar) #%% Setup acquisition and send pattern to the DMD -setup_version = 'setup_v1.3' -data_folder_name = '2022-09-02_optical_tuning' -data_name = 'DMD_12cm_f75mm_6.5cm_f50mm_6.5cm_f30mm_1.5cm_x20_0mm_diam1.5mm' +setup_version = 'setup_v1.3.1' +Np = 32 # Number of pixels in one dimension of the image (image: NpxNp) +zoom = 1 # Numerical zoom applied in the DMD +ti = 2 # Integration time of the spectrometer +scan_mode = 'Walsh_inv' #'Raster_inv' #'Walsh' #'Raster' # +data_folder_name = '2023-05-03_diffraction_unlimited' +data_name = 'cat_' + scan_mode + '_im_'+str(Np)+'x'+str(Np)+'_ti_'+str(ti)+'ms_zoom_x'+str(zoom) -camPar.acq_mode = 'snapshot'#'video' # +camPar.acq_mode = 'video' #'snapshot'# camPar.vidFormat = 'avi' #'bin'# camPar.insert_patterns = 0 # 0: no insertion / 1: insert white patterns for the camera -camPar.gate_period = 16 # a multiple of the integration time of the spectro, between [1 - 16] +camPar.gate_period = 1 # a multiple of the integration time of the spectro, between [2 - 16] (2: insert one white pattern between each pattern) camPar.black_pattern_num = 1 # insert the picture number (in the pattern_source folder) of the pattern you want to insert - all_path = func_path(data_folder_name, data_name) metadata = MetaData( output_directory = all_path.subfolder_path, - pattern_order_source = 'C:/openspyrit/spas/stats/pattern_order.npz',#'../communication/raster.txt',# - pattern_source = 'C:/openspyrit/spas/Patterns/PosNeg/DMD_Walsh_64x64',#'../Patterns/RasterScan_64x64',# - pattern_prefix = 'Walsh_64x64',#'RasterScan_64x64_1',# + pattern_order_source = 'C:/openspyrit/spas/stats/pattern_order_' + scan_mode + '_' + str(Np) + 'x' + str(Np) + '.npz', + pattern_source = 'C:/openspyrit/spas/Patterns/Zoom_x' + str(zoom) + '/' + scan_mode + '_' + str(Np) + 'x' + str(Np), + pattern_prefix = scan_mode + '_' + str(Np) + 'x' + str(Np), experiment_name = data_name, - light_source = 'White LED light',#'HgAr multilines Source (HG-1 Oceanoptics)',#'Nothing',# - object = 'The cat',#'Nothing',#'USAF',#'Star Sector',#'Nothing' - filter = 'Diffuser + HighPass_500nm + LowPass_750nm',#' linear colored filter + OD#0',#'Nothing',# - description = 'Test') + light_source = 'Zeiss KL2500 white lamp',#'White LED light',#'LED Laser 385nm + optical fiber 600µm, P = 30 mW',#'the sun',#'IKEA lamp 10W LED1734G10',#or BlueLaser 50 mW',#' (74) + 'Bioblock power: II',#'HgAr multilines Source (HG-1 Oceanoptics)',#'Nothing',# + object = 'Cat',#two little tubes containing PpIX at 634 and 620 state',#'Apple',#',#'USAF',#'Nothing''color checker' + filter = 'BandPass filter 560nm Dl=10nm',#'HighPass_500nm + LowPass_750nm',# + optical density = 0.1', #'Nothing',#'Diffuser + HighPass_500nm + LowPass_750nm',##'Microsope objective x40',#'' linear colored filter + OD#0',#'Nothing',# + description = 'numerical pinhole') acquisition_parameters = AcquisitionParameters( pattern_compression = 1.0, - pattern_dimension_x = 64, - pattern_dimension_y = 64) + pattern_dimension_x = Np, + pattern_dimension_y = Np) spectrometer_params, DMD_params, camPar = setup_2arms( spectrometer = spectrometer, @@ -76,36 +81,9 @@ metadata = metadata, acquisition_params = acquisition_parameters, DMD_output_synch_pulse_delay = 42, - integration_time = 1) -#%% Setup reconstruction -network_params = ReconstructionParameters( - img_size = 64, - CR = 1024, - denoise = True, - epochs = 40, - learning_rate = 1e-3, - step_size = 20, - gamma = 0.2, - batch_size = 256, - regularization = 1e-7, - N0 = 50.0, - sig = 0.0, - arch_name = 'c0mp') - -cov_path = 'C:/openspyrit/spas/stats/new-nicolas/Cov_64x64.npy' -mean_path = 'C:/openspyrit/spas/stats/new-nicolas/Average_64x64.npy' -model_root = 'C:/openspyrit/spas/models/new-nicolas/' -H = wh.walsh2_matrix(64)/64 - -model, device = setup_reconstruction(cov_path, mean_path, H, model_root, network_params) -noise = load_noise('C:/openspyrit/spas/noise-calibration/fit_model2.npz') - -reconstruction_params = { - 'model' : model, - 'device' : device, - 'batches': 1, - 'noise' : noise} + integration_time = ti) #%% Acquire +time.sleep(0) if camPar.acq_mode == 'video': spectral_data = acquire_2arms( ava = spectrometer, @@ -128,16 +106,43 @@ acquisition_params = acquisition_parameters, repetitions = 1, reconstruct = False) -#%% Reconstruction without NN -Q = wh.walsh2_matrix(64) -GT = reconstruction_hadamard(acquisition_parameters.patterns, 'walsh', Q, spectral_data) -plot_reco_without_NN(acquisition_parameters, GT, Q, all_path.had_reco_path, - all_path.fig_had_reco_path) -#%% Reconstruct with NN -plot_reco_with_NN(acquisition_parameters, network_params, spectral_data, noise, - model, device) +#%% Hadamard Reconstruction +Q = wh.walsh2_matrix(Np) +GT = reconstruction_hadamard(acquisition_parameters.patterns, 'walsh', Q, spectral_data, Np) +plot_reco_without_NN(acquisition_parameters, GT, Q, all_path) +#%% Neural Network Reconstruction +t0 = time.time() +network_param = ReconstructionParameters( + # Reconstruction network + M = 64*64, # Number of measurements + img_size = 128, # Image size + arch = 'dc-net', # Main architecture + denoi = 'unet', # Image domain denoiser + subs = 'rect', # Subsampling scheme + + # Training + data = 'imagenet', # Training database + N0 = 10, # Intensity (max of ph./pixel) + + # Optimisation (from train2.py) + num_epochs = 30, # Number of training epochs + learning_rate = 0.001, # Learning Rate + step_size = 10, # Scheduler Step Size + gamma = 0.5, # Scheduler Decrease Rate + batch_size = 256, # Size of the training batch + regularization = 1e-7 # Regularisation Parameter + ) + +cov_path = 'C:/openspyrit/stat/ILSVRC2012_v10102019/Cov_8_128x128.npy' +model_folder = 'C:/openspyrit/models/' +model, device = setup_reconstruction(cov_path, model_folder, network_param) +meas = reorder_subsample(spectral_data.T, acquisition_parameters, network_param) # Reorder and subsample +reco = reconstruct(model, device, meas) # Reconstruction +plot_reco_with_NN(acquisition_parameters, reco, all_path) +print('elapsed time = ' + str(round(time.time()-t0)) + ' s') #%% transfer data to girder transfer_data_2arms(metadata, acquisition_parameters, spectrometer_params, DMD_params, camPar, - setup_version, data_folder_name, data_name) + setup_version, data_folder_name, data_name, upload_metadata = 1) #%% Disconnect -disconnect_2arms(spectrometer, DMD, camPar) \ No newline at end of file +disconnect_2arms(spectrometer, DMD, camPar) + diff --git a/scripts/reconstruct_and_transfer_data_to_girder.py b/scripts/reconstruct_and_transfer_data_to_girder.py new file mode 100644 index 0000000..38020dc --- /dev/null +++ b/scripts/reconstruct_and_transfer_data_to_girder.py @@ -0,0 +1,577 @@ +# -*- coding: utf-8 -*- +__author__ = 'Guilherme Beneti Martins' + + +#%% Package +import os +import numpy as np +import csv +import time +import pickle +import girder_client + +import spyrit.misc.walsh_hadamard as wh +from spas.visualization import plot_reco_without_NN, plot_reco_with_NN +from spas.metadata import func_path, read_metadata +from spas.transfer_data_to_girder import transfer_data_2arms, transfer_data +from spas import reconstruct +from spas.reconstruction import reconstruction_hadamard +from spas.reconstruction_nn import ReconstructionParameters, setup_reconstruction, reorder_subsample + + + +#%% INPUT +delete_old_fig_init = 1 +read_data_init = 1 +had_reco_init = 1 +nn_reco_init = 1 +check_data_exist_in_girder_init = 1 +tranfer_init = 1 +#%%############################# init ############################ +inc = 0 +inc_error = 0 +list_error = [] +t_tot_0 = time.time() +########################## BEGINNING ############################ +# setup_version = 'setup_v1.3.1' +#data_folder_name = '2023-03-17_test_Leica_microscope_HCL_Bron' +data_folder_name_list = os.listdir('../data/') +# data_folder_name_list = ['2021-04-27-reconstruction']# ['2021-06-01-spectral-resolution2']#'2023-03-28_test_transfer_image'#'2023-03-17_test_flip_image'# +# data_folder_name_list = ['2022-09-19_Thorlabs_box'] +for data_folder_name in data_folder_name_list: + data_file_list = os.listdir('../data/' + data_folder_name) + # data_file_list = ['red_and_black_ink_im_64x64_ti_20ms_zoom_x2'] + # data_file_list = ['Thorlabs_box_ti_10ms_without_telecentric']#['PpIX_634_laser_385nm_im_32x32_ti_200ms_zoom_x1']# + # data_file_list = ['reconstruction-diffuser'] + ############################ CSV file ########################## + csv_file_path = 'data_in_Girder/data.csv' + csv_exist = os.path.isfile(csv_file_path) + fieldnames = ['setup_version', 'data_folder_name', 'data_name', 'transfered_to_girder', 'had_reco', 'nn_reco', 'check_data_exist_in_girder', 'delete_old_fig'] + + for data_name in data_file_list: + if data_name != 'acq-blue' and data_name != 'acq-gray' and data_name != 'acq-red': + ################## possibly to change initial input ################## + delete_old_fig = delete_old_fig_init + read_data = read_data_init + had_reco = had_reco_init + nn_reco = nn_reco_init + check_data_exist_in_girder = check_data_exist_in_girder_init + tranfer = tranfer_init + ############################ path ################################### + print('data folder : ' + data_folder_name) + print(' -' + data_name) + all_path = func_path(data_folder_name, data_name) + ###### check if data has been already traited via the CSV file ###### + duplicate = 0 + if csv_exist == True: + with open(csv_file_path, newline='') as csvfile: + reader = csv.DictReader(csvfile) + for row in reader: + # print(row[fieldnames[0]], row[fieldnames[1]], row[fieldnames[2]], row[fieldnames[3]]) + csv_setup_version = row[fieldnames[0]] + csv_data_folder_name = row[fieldnames[1]] + csv_data_name = row[fieldnames[2]] + csv_transfered_to_girder = row[fieldnames[3]] + csv_had_reco = row[fieldnames[4]] + csv_nn_reco = row[fieldnames[5]] + csv_check_data_exist_in_girder = row[fieldnames[6]] + csv_delete_old_fig = row[fieldnames[7]] + if csv_data_folder_name == data_folder_name and csv_data_name == data_name: + duplicate = 1 + if csv_delete_old_fig == '1': + delete_old_fig = 0 + if csv_had_reco == '1': + read_data = 0 + had_reco = 0 + if csv_nn_reco == '1': + nn_reco = 0 + if csv_check_data_exist_in_girder == '1': + check_data_exist_in_girder = 0 + if csv_transfered_to_girder == '1': + tranfer = 0 + ###################### delete old figures ########################### + if delete_old_fig == 1: + fig_list = os.listdir(all_path.overview_path) + print('---------- delete old figures ----------') + for fig in fig_list: + if fig.find('HAD_RECO') >=0 or fig.find('NN_RECO') >=0: + print(fig) + os.remove(all_path.overview_path + '/' + fig) + ########################## read raw data ########################### + if read_data == 1: + print('---------- read data ----------') + file = np.load(all_path.data_path+'_spectraldata.npz') + try: + spectral_data = file['spectral_data'] + # print('npz item : spectral_data') + # for k in file.data_name: + # print(k) + except: + spectral_data = file['arr_0'] + print('npz item : arr_0') + + + ########################### read metadata ########################### + metadata_path = all_path.data_path + '_metadata.json' + metadata, acquisition_parameters, spectrometer_parameters, DMD_parameters = read_metadata(metadata_path) + + if (acquisition_parameters.pattern_amount == 10000 or data_name.find('raster') >= 0 or data_name.find('Raster') >= 0 or + data_name.find('test_1_ti_4_tc_0.09_2') >= 0 or data_folder_name.find('Kalman') >= 0 or data_folder_name.find('Antonio') >= 0 or + data_folder_name.find('diffraction_unlimited') >= 0 or data_folder_name.find('128x128') >= 0 or data_name.find('128x128') >= 0 or + data_name.find('32x32') >= 0 or data_name.find('16x16') >= 0 or data_name.find('8x8') >= 0): + had_reco = 0 + nn_reco = 0 + tranfer = 0 + + wavelengths = acquisition_parameters.wavelengths + Np = acquisition_parameters.pattern_dimension_x + else: + had_reco = 0 + nn_reco = 0 + check_data_exist_in_girder = 0 + #################### Hadamard reconstruction ####################### + if had_reco == 1: + print('---------- had reco ----------') + # subsampling + nsub = 1 + meas_size = acquisition_parameters.pattern_dimension_x * acquisition_parameters.pattern_dimension_y * 2 + M_sub = spectral_data[:meas_size//nsub,:] + patterns_sub = acquisition_parameters.patterns[:meas_size//nsub] + + ### Hadamard Reconstruction + Q = wh.walsh2_matrix(Np) + GT = reconstruction_hadamard(patterns_sub, 'walsh', Q, M_sub, Np) + # GT = np.fliplr(GT) + plot_reco_without_NN(acquisition_parameters, GT, Q, all_path) + + ######################### read cam metadata ######################## + try: + cam_metadata_path = all_path.data_path + '_metadata_cam.pkl' + + file = open(cam_metadata_path,'rb') + cam_metadata = pickle.load(file) + file.close() + metadata_cam = 1 + except: + print('metada of the cam does not exist') + cam_metadata = [] + metadata_cam = 0 + + camPar = cam_metadata + + #%% Neural Network Reconstruction + if nn_reco == 1: + print('---------- nn reco ----------') + t0 = time.time() + network_param = ReconstructionParameters( + # Reconstruction network + M = round(meas_size/2), #64*64, # Number of measurements + img_size = 128, # Image size + arch = 'dc-net', # Main architecture + denoi = 'unet', # Image domain denoiser + subs = 'rect', # Subsampling scheme + + # Training + data = 'imagenet', # Training database + N0 = 10, # Intensity (max of ph./pixel) + + # Optimisation (from train2.py) + num_epochs = 30, # Number of training epochs + learning_rate = 0.001, # Learning Rate + step_size = 10, # Scheduler Step Size + gamma = 0.5, # Scheduler Decrease Rate + batch_size = 256, # Size of the training batch + regularization = 1e-7 # Regularisation Parameter + ) + + cov_path = 'C:/openspyrit/stat/ILSVRC2012_v10102019/Cov_8_128x128.npy' + model_folder = 'C:/openspyrit/models/' + model, device = setup_reconstruction(cov_path, model_folder, network_param) + # meas = reorder_subsample(spectral_data.T, acquisition_parameters, network_param) # Reorder and subsample + # reco = reconstruct(model, device, meas) # Reconstruction + plot_reco_with_NN(acquisition_parameters, spectral_data, model, device, network_param, all_path) + print('elapsed time = ' + str(round(time.time()-t0)) + ' s') + #%% check if data exist in the Pilot warehosue + transfer_matched = 0 + upload_metadata = 0 + if check_data_exist_in_girder == 1: + print('---------- check if data exist in Girder ----------') + ########################### Girder info ################################# + url = 'https://pilot-warehouse.creatis.insa-lyon.fr/api/v1' + collectionId = '6140ba6929e3fc10d47dbe3e'# collection_name = 'spc' + txt_file = open('C:/private/no_name.txt', 'r', encoding='utf8') + apiKey = txt_file.read() + txt_file.close() + #################### Girder authentification ############################ + gc = girder_client.GirderClient(apiUrl=url) # Generate the warehouse client + gc.authenticate(apiKey=apiKey) # Authentication to the warehouse + ####################### scan folders inside Girder ##################### + girder_data_folder_id = '6149c3ce29e3fc10d47dbffb' + version_list = gc.listFolder(girder_data_folder_id, 'folder') + for version_folder in version_list: + version_folder_id = version_folder['_id'] + data_folder_list = gc.listFolder(version_folder_id, 'folder') + for girder_data_folder in data_folder_list: + girder_data_folder_name = girder_data_folder['name'] + girder_data_folder_id = girder_data_folder['_id'] + girder_data_list = gc.listFolder(girder_data_folder_id, 'folder') + for girder_data_file in girder_data_list: + girder_data_file_name = girder_data_file['name'] + girder_data_file_id = girder_data_file['_id'] + if (girder_data_folder_name == data_folder_name) and (girder_data_file_name == data_name): + transfer_matched = 1 + setup_version = version_folder['name'] + print('data in Girder and in local matched') + print('------------ Setup Version : ' + setup_version + ' --------------') + print(' * ' + data_folder_name) + print(' - ' + data_name) + if data_folder_name == '2023-03-17_test_Leica_microscope_HCL_Bron' or data_folder_name == '2023-04-05_PpIX_at_lab_to_compare_to_hospital' or data_folder_name == '2023-04-07_PpIX_at_lab_to_compare_to_hospital': + upload_metadata = 1 + transfer_matched = 1 + else: + upload_metadata = 0 + + break + # else: + # print('no matching') + + if transfer_matched == 1: + break + if transfer_matched == 1: + break + if transfer_matched == 0: + inc_error = inc_error + 1 + list_error.append(data_folder_name + '/' + data_name) + + #%% transfer data to girder + transfer_done = 0 + if tranfer == 1 and transfer_matched == 1: + print('---------- transfer data ----------') + if metadata_cam == 0: + transfer_data(metadata, acquisition_parameters, spectrometer_parameters, DMD_parameters, + setup_version, data_folder_name, data_name, upload_metadata) + else: + transfer_data_2arms(metadata, acquisition_parameters, spectrometer_parameters, DMD_parameters, camPar, + setup_version, data_folder_name, data_name, upload_metadata) + + transfer_done = 1 + + + #%% Write dataLog in csv file + # save csv data + if duplicate == 0: + print('---------- write CSV file ----------') + try: + # if transfer_matched = 0, setup_version is not defined + setup_ver = setup_version + except: + setup_ver = '' + + rows = [ + {'setup_version': setup_ver, + 'data_folder_name': data_folder_name, + 'data_name': data_name, + 'transfered_to_girder': transfer_matched, + 'had_reco': had_reco, + 'nn_reco': nn_reco, + 'check_data_exist_in_girder': check_data_exist_in_girder, + 'delete_old_fig': delete_old_fig} + ] + + with open(csv_file_path, 'a', encoding='UTF8', newline='') as f: + writer = csv.DictWriter(f, fieldnames=fieldnames) + if csv_exist == False: + writer.writeheader() + writer.writerows(rows) + else: + print(' ==> data already traited !!') + + else: + print('file : ' + data_name + ' is not an hadamard acq') +print('error number (data do not match) = ' + str(inc_error)) +print('total elapsed time = ' + str(round(time.time()-t_tot_0)) + ' s') +# #%%####################### spectra ##################################### +# from scipy.signal import savgol_filter + +# plot_raw_spectrum = 0 +# plot_smooth_spectra = 1 +# plot_norm_smooth_spectra = 1 + + +# if data_name == 'Painting_Don_Quichotte_64x64_ti_150ms': +# y1 = np.mean(np.mean(GT[8:11,29:34,:], axis=1), axis=0) #'white ref' +# y2 = np.mean(np.mean(GT[25:28,27:30,:], axis=1), axis=0) #'base de la queue' +# y3 = np.mean(GT[22,33:37,:], axis=0) #'milieu de la queue' +# y4 = np.mean(np.mean(GT[17:21,19:24,:], axis=1), axis=0) #'background' +# elif data_name == 'Painting_St-Tropez_64x64_ti_100ms': +# y1 = np.mean(np.mean(GT[19:33,54:57,:], axis=1), axis=0) #'white ref' +# y2 = np.mean(np.mean(GT[15:19,21:28,:], axis=1), axis=0) #'left house' +# y3 = np.mean(np.mean(GT[20:22,41:50,:], axis=1), axis=0) #'right house' +# y4 = np.mean(np.mean(GT[45:56,15:33,:], axis=1), axis=0) #'boat' + + +# window_size = 51 +# polynomial_order = 4 +# ysm1 = savgol_filter(y1, window_size, polynomial_order) +# ysm2 = savgol_filter(y2, window_size, polynomial_order) +# ysm3 = savgol_filter(y3, window_size, polynomial_order) +# ysm4 = savgol_filter(y4, window_size, polynomial_order) + +# if plot_raw_spectrum == 1: +# plt.figure() +# plt.plot(wavelengths, y1, color = 'blue') +# plt.plot(wavelengths, ysm1, color = 'red') +# plt.title('white ref') +# plt.grid() + +# plt.figure() +# plt.plot(wavelengths, y2, color = 'blue') +# plt.plot(wavelengths, ysm2, color = 'red') +# plt.title('base de la queue') +# plt.grid() + +# plt.figure() +# plt.plot(wavelengths, y3, color = 'blue') +# plt.plot(wavelengths, ysm3, color = 'red') +# plt.title('milieu de la queue') +# plt.grid() + +# plt.figure() +# plt.plot(wavelengths, y4, color = 'blue') +# plt.plot(wavelengths, ysm4, color = 'red') +# plt.title('background') +# plt.grid() + +# if plot_smooth_spectra == 1: +# plt.figure() +# plt.plot(wavelengths, ysm1, color = 'green') +# plt.plot(wavelengths, ysm2, color = 'blue') +# plt.plot(wavelengths, ysm3, color = 'red') +# plt.plot(wavelengths, ysm4, color = 'black') +# plt.legend(['1', '2', '3', '4']) +# plt.grid() + +# if plot_norm_smooth_spectra == 1: +# cut = 10 +# ym1 = ysm1[cut:-cut]/np.amax(ysm1[cut:-cut]) +# ym2 = ysm2[cut:-cut]/ym1 +# ym3 = ysm3[cut:-cut]/ym1 +# ym4 = ysm4[cut:-cut]/ym1 + +# plt.figure() +# plt.plot(wavelengths[cut:-cut], ym2, color = 'blue') +# plt.plot(wavelengths[cut:-cut], ym3, color = 'red') +# plt.plot(wavelengths[cut:-cut], ym4, color = 'black') +# plt.legend(['2', '3', '4']) +# plt.grid() + +# if data_name == 'Falcon_620_WhiteLight_OFF_BlueLaser_ON_im_32x32_Zoom_x1_ti_1000ms#_tc_10.0ms': +# # 620nm +# plt.figure() +# plt.plot(wavelengths, np.mean(np.mean(GT[3:26,5:25,:], axis=1), axis=0)) +# plt.axvline(x = 620, color = 'b', label = 'axvline - full height') +# plt.axvline(x = 634, color = 'r', label = 'axvline - full height') +# plt.title('620 nm') +# plt.grid() + +# if data_name == 'Falcon_634_WhiteLight_OFF_BlueLaser_ON_im_64x64_Zoom_x1_ti_50ms#_tc_4.619ms': +# plt.figure() +# plt.plot(wavelengths, np.mean(np.mean(GT[15:45,15:45,:], axis=1), axis=0)) +# plt.axvline(x = 620, color = 'b', label = 'axvline - full height') +# plt.axvline(x = 634, color = 'r', label = 'axvline - full height') +# plt.title('634 nm') +# plt.grid() + +# plt.figure() +# plt.plot(wavelengths, GT[32,32,:]) +# plt.axvline(x = 620, color = 'b', label = 'axvline - full height') +# plt.axvline(x = 634, color = 'r', label = 'axvline - full height') +# plt.title('634 nm - single pixel') +# plt.grid() + +# if data_name == 'WhiteLight_OFF_BlueLaser_ON_im_64x64_Zoom_x1_ti_100ms#_tc_4.619ms': +# # 620nm +# plt.figure() +# plt.plot(wavelengths, np.mean(np.mean(GT[6:22,19:32,:], axis=1), axis=0)) +# plt.axvline(x = 620, color = 'b', label = 'axvline - full height') +# plt.axvline(x = 634, color = 'r', label = 'axvline - full height') +# plt.title('620 nm') +# plt.grid() + +# # 634nm +# plt.figure() +# plt.plot(wavelengths, np.mean(np.mean(GT[20:50,37:54,:], axis=1), axis=0)) +# plt.axvline(x = 620, color = 'b', label = 'axvline - full height') +# plt.axvline(x = 634, color = 'r', label = 'axvline - full height') +# plt.title('634 nm') +# plt.grid() + +# plt.figure() +# plt.plot(wavelengths, GT[25,44,:]) +# plt.axvline(x = 620, color = 'b', label = 'axvline - full height') +# plt.axvline(x = 634, color = 'r', label = 'axvline - full height') +# plt.title('634 nm - single pixel') +# plt.grid() + +# # S0 +# plt.figure() +# plt.plot(wavelengths, np.mean(np.mean(GT[27:52,8:30,:], axis=1), axis=0)) +# plt.axvline(x = 620, color = 'b', label = 'axvline - full height') +# plt.axvline(x = 634, color = 'r', label = 'axvline - full height') +# plt.title('S_0') +# plt.grid() +######################### subsampling ####################################### +# ============================================================================= +# N = 64 +# nsub = 2 +# M_sub = M[:8192//nsub,:] +# acquisition_parameters.patterns_sub = acquisition_parameters.patterns[:8192//nsub] +# GT_sub = reconstruction_hadamard(acquisition_parameters.patterns_sub, 'walsh', Q, M_sub) +# F_bin_sub, wavelengths_bin, bin_width = spectral_binning(GT_sub.T, acquisition_parameters.wavelengths, 530, 730, 8) +# +# +# +# plot_color(F_bin_sub, wavelengths_bin) +# plt.savefig(fig_had_reco_path + '_wavelength_binning_subsamplig=' + str(nsub) + '.png') +# plt.show() +# ============================================================================= + +# ============================================================================= +# plt.figure +# plt.imshow(GT[:,:,0]) +# plt.title(f'lambda = {wavelengths[0]:.2f} nm') +# plt.savefig(fig_had_reco_path + '_' + f'lambda = {wavelengths[0]:.2f} nm.png') +# plt.show() +# +# plt.figure +# plt.imshow(GT[:,:,410]) +# plt.title(f'lambda = {wavelengths[410]:.2f} nm') +# plt.savefig(fig_had_reco_path + '_' + f'lambda = {wavelengths[410]:.2f} nm.png') +# plt.show() +# +# plt.figure +# plt.imshow(GT[:,:,820]) +# plt.title(f'lambda = {wavelengths[820]:.2f} nm') +# plt.savefig(fig_had_reco_path + '_' + f'lambda = {wavelengths[820]:.2f} nm.png') +# plt.show() +# +# plt.figure +# plt.imshow(GT[:,:,1230]) +# plt.title(f'lambda = {wavelengths[1230]:.2f} nm') +# plt.savefig(fig_had_reco_path + '_' + f'lambda = {wavelengths[1230]:.2f} nm.png') +# plt.show() +# +# plt.figure +# plt.imshow(GT[:,:,2047]) +# plt.title(f'lambda = {wavelengths[2047]:.2f} nm') +# plt.savefig(fig_had_reco_path + '_' + f'lambda = {wavelengths[2047]:.2f} nm.png') +# plt.show() +# +# plt.figure +# plt.imshow(np.sum(GT,axis=2)) +# plt.title('Sum of all wavelengths') +# plt.savefig(fig_had_reco_path + '_sum_of_wavelengths.png') +# plt.show() +# +# plt.figure +# plt.scatter(wavelengths, np.mean(np.mean(GT,axis=1),axis=0)) +# plt.grid() +# plt.xlabel('Lambda (nm)') +# plt.title('Spectral view in the spatial mean') +# plt.savefig(fig_had_reco_path + '_spectral_axe_of_the_hypercube.png') +# plt.show() +# +# indx = np.where(GT == np.max(GT)) +# sp = GT[indx[0],indx[1],:] +# plt.figure +# plt.scatter(wavelengths, sp.T) +# plt.grid() +# plt.xlabel('Lambda (nm)') +# plt.title('Spectral view of the max intensity') +# plt.savefig(fig_had_reco_path + '_spectral_axe_of_max_intensity.png') +# plt.show() +# ============================================================================= +# #%% Reconstruct with NN + +# #%% Setup reconstruction +# network_params = ReconstructionParameters( +# img_size = Np, +# CR = 1024, +# denoise = True, +# epochs = 40, +# learning_rate = 1e-3, +# step_size = 20, +# gamma = 0.2, +# batch_size = 256, +# regularization = 1e-7, +# N0 = 50.0, +# sig = 0.0, +# arch_name = 'c0mp') + +# cov_path = 'C:/openspyrit/spas/stats/Cov_'+str(Np)+'x'+str(Np)+'.npy' +# mean_path = 'C:/openspyrit/spas/stats/Average_'+str(Np)+'x'+str(Np)+'.npy' +# model_root = 'C:/openspyrit/spas/models/new-nicolas/' +# H = wh.walsh2_matrix(Np)/Np + +# model, device = setup_reconstruction(cov_path, mean_path, H, model_root, network_params) +# noise = load_noise('C:/openspyrit/spas/noise-calibration/fit_model2.npz') + +# reconstruction_params = { +# 'model' : model, +# 'device' : device, +# 'batches': 1, +# 'noise' : noise} + +# # network_params = ReconstructionParameters( +# # img_size=64, +# # CR=1024, +# # denoise=True, +# # epochs=40, +# # learning_rate=1e-3, +# # step_size=20, +# # gamma=0.2, +# # batch_size=256, +# # regularization=1e-7, +# # N0=50.0, +# # sig=0.0, +# # arch_name='c0mp',) + +# # cov_path = '../stats/new-nicolas/Cov_64x64.npy' +# # mean_path = '../stats/new-nicolas/Average_64x64.npy' +# # H_path = '../stats/new-nicolas/H.npy' +# # model_root = '../models/new-nicolas/' + +# # model, device = setup_reconstruction(cov_path, mean_path, H_path, model_root, network_params) +# # noise = load_noise('../noise-calibration/fit_model2.npz') + +# # reconstruction_params = { +# # 'model': model, +# # 'device': device, +# # 'batches': 1, +# # 'noise': noise, +# # } + +# F_bin, wavelengths_bin, bin_width, noise_bin = spectral_binning(M.T, acquisition_parameters.wavelengths, 530, 730, 8, 0, noise) +# recon = reconstruct(model, device, F_bin[0:8192//4,:], 1, noise_bin) +# plot_color(recon, wavelengths_bin) +# plt.savefig(nn_reco_path + '_reco_wavelength_binning.png') +# plt.show() + + +# #%% transfer data to girder +# transf.transfer_data_to_girder(metadata, acquisition_parameters, spectrometer_params, DMD_params, setup_version, data_folder_name, data_name) +#%% delete plots +# Question = input("Do you want to delete the figures yes [y] ? ") +# if Question == ("y") or Question == ("yes"): +# shutil.rmtree(overview_path) +# print ("==> figures deleted") + + + + + + + + + + diff --git a/spas/acquisition.py b/spas/acquisition.py index 5fb6eaa..7701c0b 100644 --- a/spas/acquisition.py +++ b/spas/acquisition.py @@ -52,16 +52,23 @@ from collections import namedtuple from pathlib import Path from multiprocessing import Process, Queue -import shutil - +import shutil import numpy as np from PIL import Image -from msl.equipment import EquipmentRecord, ConnectionRecord, Backend -from msl.equipment.resources.avantes import MeasureCallback, Avantes -from ALP4 import ALP4, ALP_FIRSTFRAME, ALP_LASTFRAME -from ALP4 import ALP_AVAIL_MEMORY, ALP_DEV_DYN_SYNCH_OUT1_GATE, tAlpDynSynchOutGate -from tqdm import tqdm +try: + from ALP4 import ALP4, ALP_FIRSTFRAME, ALP_LASTFRAME + from ALP4 import ALP_AVAIL_MEMORY, ALP_DEV_DYN_SYNCH_OUT1_GATE, tAlpDynSynchOutGate +except: + class ALP4: + pass +try: + from msl.equipment import EquipmentRecord, ConnectionRecord, Backend + from msl.equipment.resources.avantes import MeasureCallback, Avantes +except: + pass + +from tqdm import tqdm from .metadata import DMDParameters, MetaData, AcquisitionParameters from .metadata import SpectrometerParameters, save_metadata, CAM, save_metadata_2arms from .reconstruction_nn import reconstruct_process, plot_recon, ReconstructionParameters @@ -638,6 +645,16 @@ def _setup_timings(DMD: ALP4, DMD_params: DMDParameters, picture_time: int, triggerInDelay=trigger_in_delay) DMD_params.update_sequence_parameters(add_illumination_time, DMD=DMD) + + +# class mytAlpDynSynchOutGate(ct.Structure): +# # For ControlType ALP_DEV_DYN_TRIG_OUT[1..3]_GATE of function AlpDevControlEx +# # Configure compiler to not insert padding bytes! (e.g. #pragma pack) +# _pack_ = 1 +# _fields_ = [("Period", ct.c_ubyte), # Period=1..16 enables output; 0: tri-state +# ("Polarity", ct.c_ubyte), # 0: active pulse is low, 1: high +# ("Gate", ct.c_ubyte * 16), +# ("byref", ct.c_ubyte * 18)] def setup(spectrometer: Avantes, DMD: ALP4, diff --git a/spas/metadata.py b/spas/metadata.py index ad38de1..97ea932 100644 --- a/spas/metadata.py +++ b/spas/metadata.py @@ -8,20 +8,33 @@ improve readability. """ -import json -import ALP4 +import json from datetime import datetime from enum import IntEnum from dataclasses import dataclass, InitVar, field from typing import Optional, Union, List, Tuple, Optional from pathlib import Path import os - -from msl.equipment.resources.avantes import MeasConfigType from dataclasses_json import dataclass_json import numpy as np import ctypes as ct import pickle +##### DLL for the DMD +try: + from ALP4 import ALP4 +except: # in the cas the DLL of the DMD is not installed + class ALP4: + pass + setattr(ALP4, 'ALP4', None) + print('DLL of the DMD not installed') +##### DLL for the spectrometer Avantes +try: + from msl.equipment.resources.avantes import MeasConfigType +except: # in the cas the DLL of the spectrometer is not installed + class MeasConfigType: + pass + MeasConfigType = None + print('DLL of the spectrometer not installed !!!') class DMDTypes(IntEnum): diff --git a/spas/reconstruction_nn.py b/spas/reconstruction_nn.py index b5d3242..3863689 100644 --- a/spas/reconstruction_nn.py +++ b/spas/reconstruction_nn.py @@ -9,45 +9,25 @@ """ from dataclasses import InitVar, dataclass, field -from enum import IntEnum from time import perf_counter_ns, sleep from typing import Tuple, Union from multiprocessing import Queue -from pathlib import Path - import math import torch import numpy as np from matplotlib import pyplot as plt -# from spyrit.core.reconstruction import Pinv_Net, DC2_Net -from spyrit.core.recon import PinvNet, DCNet # modified by LMW 30/03/2023 - -# from spyrit.core.training import load_net -from spyrit.core.train import load_net # modified by LMW 30/03/2023 - +from spyrit.core.recon import PinvNet, DCNet +from spyrit.core.train import load_net from spyrit.misc.statistics import Cov2Var from spyrit.misc.walsh_hadamard import walsh2_matrix - -# from spyrit.core.Acquisition import Acquisition_Poisson_approx_Gauss -from spyrit.core.noise import Poisson # modified by LMW 30/03/2023 - -# from spyrit.core.Forward_Operator import Forward_operator_Split_ft_had -from spyrit.core.meas import HadamSplit # modified by LMW 30/03/2023 - -# from spyrit.core.Preprocess import Preprocess_Split_diag_poisson -from spyrit.core.prep import SplitPoisson # modified by LMW 30/03/2023 -# il y a aussi la classe SplitRowPoisson ??? - -# from spyrit.core.Data_Consistency import Generalized_Orthogonal_Tikhonov #, Pinv_orthogonal -from spyrit.core.recon import TikhonovMeasurementPriorDiag # modified by LMW 30/03/2023 - -# from spyrit.core.neural_network import Unet #, Identity -from spyrit.core.nnet import Unet # modified by LMW 30/03/2023 - +from spyrit.core.noise import Poisson +from spyrit.core.meas import HadamSplit +from spyrit.core.prep import SplitPoisson +from spyrit.core.recon import TikhonovMeasurementPriorDiag +from spyrit.core.nnet import Unet from spyrit.misc.sampling import Permutation_Matrix, reorder - from spas.noise import noiseClass from spas.metadata import AcquisitionParameters @@ -184,7 +164,7 @@ def setup_reconstruction(cov_path: str, net_order = 'var' bs = 256 - # net_suffix = f'N0_{network_params.N0}_N_{network_params.img_size}_M_{network_params.M}_epo_30_lr_0.001_sss_10_sdr_0.5_bs_{bs}_reg_1e-07_light' + # net_suffix = f'N0_{network_params.N0}_N_{network_params.img_size}_M_{network_params.M}_epo_30_lr_0.001_sss_10_sdr_0.5_bs_{bs}_reg_1e-07' net_suffix = f'N0_{network_params.N0}_N_{network_params.img_size}_M_{network_params.M}_epo_30_lr_0.001_sss_10_sdr_0.5_bs_{bs}_reg_1e-07_light' net_folder= f'{net_arch}_{net_denoi}_{net_data}/' @@ -345,7 +325,6 @@ def reconstruct(model: Union[PinvNet, DCNet], proportion = spectral_data.shape[0]//batches # Amount of wavelengths per batch # img_size = model.Acq.FO.h # image assumed to be square - # img_size = 128 # Modified by LMW 30/03/2023 img_size = model.Acq.meas_op.h recon = np.zeros((spectral_data.shape[0], img_size, img_size)) diff --git a/spas/transfer_data_to_girder.py b/spas/transfer_data_to_girder.py index 456d54b..3e9f53d 100644 --- a/spas/transfer_data_to_girder.py +++ b/spas/transfer_data_to_girder.py @@ -139,10 +139,16 @@ def transfer_data_2arms(metadata, acquisition_parameters, spectrometer_params, D setup_version, data_folder_name, data_name, upload_metadata): #unwrap structure into camPar - camPar.AOI_X = camPar.rectAOI.s32X.value - camPar.AOI_Y = camPar.rectAOI.s32Y.value - camPar.AOI_Width = camPar.rectAOI.s32Width.value - camPar.AOI_Height = camPar.rectAOI.s32Height.value + try: + camPar.AOI_X = camPar.rectAOI.s32X.value + camPar.AOI_Y = camPar.rectAOI.s32Y.value + camPar.AOI_Width = camPar.rectAOI.s32Width.value + camPar.AOI_Height = camPar.rectAOI.s32Height.value + except: + camPar['AOI_X'] = camPar['rectAOI'].s32X.value + camPar['AOI_Y'] = camPar['rectAOI'].s32Y.value + camPar['AOI_Width'] = camPar['rectAOI'].s32Width.value + camPar['AOI_Height'] = camPar['rectAOI'].s32Height.value #%%########################## Girder info ################################# url = 'https://pilot-warehouse.creatis.insa-lyon.fr/api/v1' collectionId = '6140ba6929e3fc10d47dbe3e'# collection_name = 'spc' diff --git a/spas/visualization.py b/spas/visualization.py index adc0456..beebe3b 100644 --- a/spas/visualization.py +++ b/spas/visualization.py @@ -6,14 +6,14 @@ from typing import Tuple, Optional import os -from spas import * -#from spas.reconstruction_nn import reconstruct +# from spas import * import numpy as np from matplotlib.colors import ListedColormap from matplotlib import pyplot as plt from mpl_toolkits.axes_grid1 import make_axes_locatable from spas.plot_spec_to_rgb_image import plot_spec_to_rgb_image -from .noise import noiseClass +from spas.noise import noiseClass +from spas.reconstruction_nn import reorder_subsample, reconstruct # Libraries for the IDS CAMERA try: @@ -483,8 +483,9 @@ def plot_reco_without_NN(acquisition_parameters, GT, Q, all_path): if not os.path.exists(had_reco_path): np.savez_compressed(had_reco_path, GT) - GT = np.rot90(GT, 1) - GT = np.rot90(GT, 1) + GT = np.rot90(GT, 2) + size_x = GT.shape[0] + size_y = GT.shape[1] F_bin, wavelengths_bin, bin_width = spectral_binning(GT.T, acquisition_parameters.wavelengths, 530, 730, 8) F_bin_rot = np.rot90(F_bin, axes=(1,2)) @@ -506,7 +507,7 @@ def plot_reco_without_NN(acquisition_parameters, GT, Q, all_path): ############### spatial view, wavelength sum ############# #plt.figure() - plt.imshow(np.sum(GT[:,:,193:877], axis=2))#[:,:,193:877] #(540-625 nm) + plt.imshow(np.mean(GT[:,:,100:-100], axis=2))#[:,:,193:877] #(540-625 nm) plt.title('Sum of all wavelengths') plt.savefig(fig_had_reco_path + '_GRAY_IMAGE_had_reco.png') plt.show() @@ -522,8 +523,8 @@ def plot_reco_without_NN(acquisition_parameters, GT, Q, all_path): plt.savefig(fig_had_reco_path + '_RGB_IMAGE_had_reco.png') plt.show() ####################### spectral view ################### - GT50 = GT[16:48,16:48,:] - GT25 = GT[24:40,24:40,:] + GT50 = GT[round(size_x/4):round(size_x*3/4), round(size_y/4):round(size_y*3/4), :] + GT25 = GT[round(size_x*3/8):round(size_x*5/8), round(size_y*3/8):round(size_y*5/8), :] plt.figure() plt.plot(acquisition_parameters.wavelengths, np.mean(np.mean(GT25,axis=1),axis=0)) plt.plot(acquisition_parameters.wavelengths, np.mean(np.mean(GT50,axis=1),axis=0)) @@ -536,40 +537,49 @@ def plot_reco_without_NN(acquisition_parameters, GT, Q, all_path): plt.show() -def plot_reco_with_NN(acquisition_parameters, reco, all_path): - - reco = reco.transpose(1, 2, 0) +def plot_reco_with_NN(acquisition_parameters, spectral_data, model, device, network_param, all_path): + + reorder_spectral_data = reorder_subsample(spectral_data.T, acquisition_parameters, network_param) + reco = reconstruct(model, device, reorder_spectral_data) # Reconstruction + reco = reco.T + reco = np.rot90(reco, 3, axes=(0,1)) nn_reco_path = all_path.nn_reco_path fig_nn_reco_path = all_path.fig_nn_reco_path if not os.path.exists(nn_reco_path): np.savez_compressed(nn_reco_path, reco) - - reco = np.rot90(reco, 1) - reco = np.rot90(reco, 1) - - F_bin, wavelengths_bin, bin_width = spectral_binning(reco.T, acquisition_parameters.wavelengths, 530, 730, 8) - F_bin_rot = np.rot90(F_bin, axes=(1,2)) - F_bin_flip = F_bin_rot[:,::-1,:] - F_bin_1px, wavelengths_bin, bin_width = spectral_slicing(reco.T, acquisition_parameters.wavelengths, 530, 730, 8) - F_bin_1px_rot = np.rot90(F_bin_1px, axes=(1,2)) - F_bin_1px_flip = F_bin_1px_rot[:,::-1,:] - ############### spatial view, wavelength bin ############# - #plt.figure() - plot_color(F_bin_flip, wavelengths_bin) - plt.savefig(fig_nn_reco_path + '_BIN_IMAGE_nn_reco.png') - plt.show() ############### spatial view, one wavelength ############# + meas_bin_1w, wavelengths_bin, _ = spectral_slicing(reorder_spectral_data, acquisition_parameters.wavelengths, 530, 730, 8) + rec = reconstruct(model, device, meas_bin_1w) # Reconstruction + rec = np.rot90(rec, 2, axes=(1,2)) + #plt.figure() - plot_color(F_bin_1px_flip, wavelengths_bin) + plot_color(rec, wavelengths_bin) plt.savefig(fig_nn_reco_path + '_SLICE_IMAGE_nn_reco.png') plt.show() - - ############### spatial view, wavelength sum ############# + + ############### spatial view, wavelength bin ############# + meas_bin, wavelengths_bin, _ = spectral_binning(reorder_spectral_data, acquisition_parameters.wavelengths, 530, 730, 8) + rec = reconstruct(model, device, meas_bin) + rec = np.rot90(rec, 2, axes=(1,2)) + #plt.figure() - plt.imshow(np.sum(reco[:,:,193:877], axis=2))#[:,:,193:877] #(540-625 nm) + plot_color(rec, wavelengths_bin) + plt.savefig(fig_nn_reco_path + '_BIN_IMAGE_nn_reco.png') + plt.show() + + ############### spatial view, wavelength sum ############# + sum_wave = np.zeros((1, reorder_spectral_data.shape[1])) + moy = np.sum(reorder_spectral_data, axis=0) + sum_wave[0, :] = moy + rec_sum = reconstruct(model, device, sum_wave) + rec_sum = rec_sum[0, :, :] + rec_sum = np.rot90(rec_sum, 2) + + # plt.figure() + plt.imshow(rec_sum)#[:,:,193:877] #(540-625 nm) plt.title('Sum of all wavelengths') plt.savefig(fig_nn_reco_path + '_GRAY_IMAGE_nn_reco.png') plt.show() @@ -577,14 +587,18 @@ def plot_reco_with_NN(acquisition_parameters, reco, all_path): ####################### RGB view ######################## print('Beging RGB convertion ...') image_arr = plot_spec_to_rgb_image(reco, acquisition_parameters.wavelengths) + image_arr = image_arr[:,::-1,:] + image_arr = np.rot90(image_arr, 2) print('RGB convertion finished') plt.figure() plt.imshow(image_arr) plt.savefig(fig_nn_reco_path + '_RGB_IMAGE_nn_reco.png') plt.show() ####################### spectral view ################### - GT50 = reco[16:48,16:48,:] - GT25 = reco[24:40,24:40,:] + size_x = reco.shape[0] + size_y = reco.shape[1] + GT50 = reco[round(size_x/4):round(size_x*3/4), round(size_y/4):round(size_y*3/4), :] + GT25 = reco[round(size_x*3/8):round(size_x*5/8), round(size_y*3/8):round(size_y*5/8), :] plt.figure() plt.plot(acquisition_parameters.wavelengths, np.mean(np.mean(GT25,axis=1),axis=0)) plt.plot(acquisition_parameters.wavelengths, np.mean(np.mean(GT50,axis=1),axis=0)) @@ -594,4 +608,5 @@ def plot_reco_with_NN(acquisition_parameters, reco, all_path): plt.legend(['25%', '50%', '100%']) plt.xlabel(r'$\lambda$ (nm)') plt.savefig(fig_nn_reco_path + '_SPECTRA_PLOT_nn_reco.png') - plt.show() \ No newline at end of file + plt.show() + From 15d95f703d665a7e3c599be5625d416f7548653f Mon Sep 17 00:00:00 2001 From: Mahieu-Williame Date: Fri, 12 May 2023 13:35:13 +0200 Subject: [PATCH 10/11] DLL import of the spectro and DMD are been placed under a TRY function to use SPAS reconstruction on Linux and MacOS plateform --- scripts/main_seq_2arms.py | 26 +-- scripts/read_video_from_monchromeIDScam.py | 209 ++++++++++++++++++ ...reconstruct_and_transfer_data_to_girder.py | 2 - spas/acquisition.py | 15 +- spas/metadata.py | 2 +- 5 files changed, 232 insertions(+), 22 deletions(-) create mode 100644 scripts/read_video_from_monchromeIDScam.py diff --git a/scripts/main_seq_2arms.py b/scripts/main_seq_2arms.py index 405736e..b7413f1 100644 --- a/scripts/main_seq_2arms.py +++ b/scripts/main_seq_2arms.py @@ -20,10 +20,10 @@ spectrometer, DMD, DMD_initial_memory, camPar = init_2arms() #%% Define the AOI # Warning, not all values are allowed for Width and Height (max: 2076x3088 | ex: 768x544) -camPar.rectAOI.s32X.value = 0 #0#1100# // X -camPar.rectAOI.s32Y.value = 0#0#800#765# // Y -camPar.rectAOI.s32Width.value = 1544 #3000#3088#3088# 768 // Width must be multiple of 8 -camPar.rectAOI.s32Height.value = 1730#2000## 1038# 544 # 2076 // Height +camPar.rectAOI.s32X.value = 1100#0 #0# // X +camPar.rectAOI.s32Y.value = 765#0#0#800# // Y +camPar.rectAOI.s32Width.value = 768# 1544 #3000#3088#3088# // Width must be multiple of 8 +camPar.rectAOI.s32Height.value = 544 #1730#2000## 1038# 2076 // Height camPar = captureVid(camPar) #%% Set Camera Parameters @@ -34,7 +34,7 @@ Gain = 0, # Gain boundary : [0 100] gain_boost = 'OFF', # set "ON" to activate gain boost, "OFF" to deactivate nGamma = 1, # Gamma boundary : [1 - 2.2] - ExposureTime = 1,# Exposure Time (ms) boudary : [0.013 - 56.221] + ExposureTime = 0.04,# Exposure Time (ms) boudary : [0.013 - 56.221] black_level = 5) # Black Level boundary : [0 255] snapshotVisu(camPar) @@ -42,14 +42,14 @@ displayVid(camPar) #%% Setup acquisition and send pattern to the DMD setup_version = 'setup_v1.3.1' -Np = 32 # Number of pixels in one dimension of the image (image: NpxNp) +Np = 64 # Number of pixels in one dimension of the image (image: NpxNp) zoom = 1 # Numerical zoom applied in the DMD -ti = 2 # Integration time of the spectrometer -scan_mode = 'Walsh_inv' #'Raster_inv' #'Walsh' #'Raster' # -data_folder_name = '2023-05-03_diffraction_unlimited' +ti = 1 # Integration time of the spectrometer +scan_mode = 'Walsh' #'Raster' #'Walsh_inv' #'Raster_inv' # +data_folder_name = '2023-05-12_test_ALP4' data_name = 'cat_' + scan_mode + '_im_'+str(Np)+'x'+str(Np)+'_ti_'+str(ti)+'ms_zoom_x'+str(zoom) -camPar.acq_mode = 'video' #'snapshot'# +camPar.acq_mode = 'snapshot'#'video' # camPar.vidFormat = 'avi' #'bin'# camPar.insert_patterns = 0 # 0: no insertion / 1: insert white patterns for the camera camPar.gate_period = 1 # a multiple of the integration time of the spectro, between [2 - 16] (2: insert one white pattern between each pattern) @@ -63,10 +63,10 @@ pattern_prefix = scan_mode + '_' + str(Np) + 'x' + str(Np), experiment_name = data_name, - light_source = 'Zeiss KL2500 white lamp',#'White LED light',#'LED Laser 385nm + optical fiber 600µm, P = 30 mW',#'the sun',#'IKEA lamp 10W LED1734G10',#or BlueLaser 50 mW',#' (74) + 'Bioblock power: II',#'HgAr multilines Source (HG-1 Oceanoptics)',#'Nothing',# + light_source = 'White LED light',#'Zeiss KL2500 white lamp',#'LED Laser 385nm + optical fiber 600µm, P = 30 mW',#'the sun',#'IKEA lamp 10W LED1734G10',#or BlueLaser 50 mW',#' (74) + 'Bioblock power: II',#'HgAr multilines Source (HG-1 Oceanoptics)',#'Nothing',# object = 'Cat',#two little tubes containing PpIX at 634 and 620 state',#'Apple',#',#'USAF',#'Nothing''color checker' - filter = 'BandPass filter 560nm Dl=10nm',#'HighPass_500nm + LowPass_750nm',# + optical density = 0.1', #'Nothing',#'Diffuser + HighPass_500nm + LowPass_750nm',##'Microsope objective x40',#'' linear colored filter + OD#0',#'Nothing',# - description = 'numerical pinhole') + filter = 'None', #'BandPass filter 560nm Dl=10nm',#'HighPass_500nm + LowPass_750nm',# + optical density = 0.1', #'Nothing',#'Diffuser + HighPass_500nm + LowPass_750nm',##'Microsope objective x40',#'' linear colored filter + OD#0',#'Nothing',# + description = 'test after changing metadata.py and acquisition.py to be used on Linux and MacOs plateform. We would like to be sure that SPAS for Windows is ok') acquisition_parameters = AcquisitionParameters( pattern_compression = 1.0, diff --git a/scripts/read_video_from_monchromeIDScam.py b/scripts/read_video_from_monchromeIDScam.py new file mode 100644 index 0000000..e571f63 --- /dev/null +++ b/scripts/read_video_from_monchromeIDScam.py @@ -0,0 +1,209 @@ +# -*- coding: utf-8 -*- +""" +Created on Tue Jan 18 17:21:21 2022 + +@author: mahieu +""" + +from matplotlib import pyplot as plt +import cv2 +import numpy as np +import pickle +from scipy.fft import fft, fftfreq + +######################## INPUT ############################### +data_folder_name = '2022-10-19_lamb_brain_motion_plus_PpIX_location_3' +data_name = 'Motion_2.5mm_WhiteLight_ON_BlueLaser_OFF_PpIX_NO_im_64x64_Zoom_x1_ti_20ms' +video_format = 'avi' +############################################################## +filename = '../data/' + data_folder_name + '/' + data_name + '/' + data_name + '_video.' + video_format +filename_out = '../data/' + data_folder_name + '/' + data_name + '/' + data_name + '_video_fps.' + video_format +if filename.find('.avi') > 1: + vid_format = 'avi' +elif filename.find('.bin') > 1: + vid_format = 'raw' +else: + vid_format = 'unknown' + print('video format unknown') +################################### read metadata cam ####### +data_path = '../data/' + data_folder_name + '/' + data_name + '/' + data_name +cam_metadata_path = data_path + '_metadata_cam.pkl' +file = open(cam_metadata_path,'rb') +cam_metadata = pickle.load(file) +file.close() + +tti = cam_metadata['time_array'] +gate_period = cam_metadata['gate_period'] +ti = tti[:len(tti)//gate_period] +delta_ti = (ti[-1]-ti[0])/len(ti) +real_fps = 1/delta_ti +#%%############################### RAW Video #################################### +if vid_format == 'raw': + ######## extract video header ################## + file = open(filename,"rb") + data_header = np.fromfile(file, dtype="uint16") + file.close() + video_header = 40 + frame_header = 48 + video_header_str = data_header[1:video_header] + + ######## extract video features ############## + Width = np.uint32(video_header_str[9]) + Height = np.uint32(video_header_str[11]) + frame_nbr = np.uint32(video_header_str[15]) + fps = float() + + ######## read frame header ############## + frame_header_str = [] + for ii in range(frame_nbr, 0, -1): + B = np.arange(start = Width*Height*(ii-1), stop = Width*Height*(ii-1)+frame_header, step = 1, dtype=np.uint32) + frame_header_str.append(B) + + ###### read real data, take off header of the video and for each frame ######## + file = open(filename,"rb") + dat = np.fromfile(file, dtype="uint8") + file.close() + ## ectract video header ## + dat3 = dat[video_header:] + ## extract frame header ## + #frame_header_str = [] + for ii in range(frame_nbr, 0, -1): + B = np.arange(start = Width*Height*(ii-1), stop = Width*Height*(ii-1)+frame_header, step = 1, dtype=np.uint32) + #frame_header_str.append(B) + dat3 = np.delete(dat3,[B]) + + ###### reshape matrix ################# + mat = np.reshape(dat3, (frame_nbr, Height, Width)) + + ## plot indiviual frames + plot_nbr = frame_nbr + if frame_nbr >= 20: + plot_nbr = 20 + + for ii in range(plot_nbr): + plt.figure() + plt.imshow(mat[ii,:,:]) + plt.colorbar(); + plt.title('frame = ' + str(ii)) + + +############################### AVI Video #################################### +elif vid_format == 'avi': + cap = cv2.VideoCapture(filename) + + Height = int(cap.get(cv2.CAP_PROP_FRAME_HEIGHT)) # float `height` + Width = int(cap.get(cv2.CAP_PROP_FRAME_WIDTH)) # float `width` + fps = cap.get(cv2.CAP_PROP_FPS) # float `fps` + frame_nbr = int(cap.get(cv2.CAP_PROP_FRAME_COUNT)) # float `total_frame_in_the_vid + # out = cv2.VideoWriter(filename_out,cv2.VideoWriter_fourcc('M','J','P','G'), real_fps, (Height,Width)) + fourcc = cv2.VideoWriter_fourcc(*'XVID') + out = cv2.VideoWriter(filename_out, fourcc, real_fps, (Width, Height)) + vid = np.zeros((Height, Width, frame_nbr)) + + for ii in range(frame_nbr): + ret, frame = cap.read() + vid[:,:,ii] = frame[:,:,0] + if ret == True: + out.write(frame) + if ii <= 15 or (ii >= frame_nbr-10 and ii <= frame_nbr-1): + plt.figure() + plt.imshow(frame[:,:,0]) + plt.colorbar(); + plt.title('frame = ' + str(ii)) + #plt.show() + + + + cap.release() + cv2.destroyAllWindows() + + + +print('') +print('height:', Height) +print('width:', Width) +print('fps:', fps) # float `fps` +print('frames count:', frame_nbr) # float `frame_count` +#%% exploit video from avi file +mean_2D = np.mean(np.mean(vid, axis = 1), axis=0) +deb = 100#0# +fin = 300#4095# + +plt.figure() +plt.plot(ti[deb:fin], mean_2D[deb:fin]) +plt.xlabel('s') +plt.title('ALL') +plt.grid() +plt.savefig(data_path + '_time_ALL.png') +##### FFT ###### +mean_2D = mean_2D - np.mean(mean_2D) +yf = fft(mean_2D) +xf = fftfreq(frame_nbr, delta_ti)[:frame_nbr//2] + +plt.figure() +plt.plot(xf, 2.0/frame_nbr * np.abs(yf[0:frame_nbr//2])) +plt.title('ALL') +plt.grid() +plt.savefig(data_path + '_fft_ALL.png') + +#%% plot part of the video +zone = np.array([[230,240,320,330],[250,300,480,520],[350,400,225,375],[230,230,320,320],[375,375,350,350]], np.int16) + +for ii in range(len(zone)): + fov = np.mean(np.mean(vid[zone[ii,:],:], axis = 1), axis=0) + + plt.figure() + plt.plot(ti[deb:fin], fov[deb:fin]) + plt.xlabel('s') + + if (zone[ii,0] == zone[ii,1]) & (zone[ii,2] == zone[ii,3]) == 1: + plt.title('pixel [' + str(zone[ii,0]) + ',' + str(zone[ii,2]) + ']') + plt.grid() + plt.savefig(data_path + '_time_pixel_[' + str(zone[ii,0]) + ',' + str(zone[ii,2]) + '].png') + else: + plt.title('zone [' + str(zone[ii,0]) + ':' + str(zone[ii,1]) + ',' + str(zone[ii,2]) + ':' + str(zone[ii,3]) + ']') + plt.grid() + plt.savefig(data_path + '_time_zone_[' + str(zone[ii,0]) + '-' + str(zone[ii,1]) + ',' + str(zone[ii,2]) + '-' + str(zone[ii,3]) + '].png') +##### FFT ###### + fov = fov - np.mean(fov) + yf = fft(fov) + xf = fftfreq(frame_nbr, delta_ti)[:frame_nbr//2] + + plt.figure() + plt.plot(xf, 2.0/frame_nbr * np.abs(yf[0:frame_nbr//2])) + if (zone[ii,0] == zone[ii,1]) & (zone[ii,2] == zone[ii,3]) == 1: + plt.title('pixel [' + str(zone[ii,0]) + ',' + str(zone[ii,2]) + ']') + plt.grid() + plt.savefig(data_path + '_fft_pixel_[' + str(zone[ii,0]) + ',' + str(zone[ii,2]) + '].png') + else: + plt.title('zone [' + str(zone[ii,0]) + ':' + str(zone[ii,1]) + ',' + str(zone[ii,2]) + ':' + str(zone[ii,3]) + ']') + plt.grid() + plt.savefig(data_path + '_fft_zone_[' + str(zone[ii,0]) + '-' + str(zone[ii,1]) + ',' + str(zone[ii,2]) + '-' + str(zone[ii,3]) + '].png') + + + + + + + + + + + + + + + + + + + + + + + + + + + + \ No newline at end of file diff --git a/scripts/reconstruct_and_transfer_data_to_girder.py b/scripts/reconstruct_and_transfer_data_to_girder.py index 38020dc..47f319f 100644 --- a/scripts/reconstruct_and_transfer_data_to_girder.py +++ b/scripts/reconstruct_and_transfer_data_to_girder.py @@ -187,8 +187,6 @@ cov_path = 'C:/openspyrit/stat/ILSVRC2012_v10102019/Cov_8_128x128.npy' model_folder = 'C:/openspyrit/models/' model, device = setup_reconstruction(cov_path, model_folder, network_param) - # meas = reorder_subsample(spectral_data.T, acquisition_parameters, network_param) # Reorder and subsample - # reco = reconstruct(model, device, meas) # Reconstruction plot_reco_with_NN(acquisition_parameters, spectral_data, model, device, network_param, all_path) print('elapsed time = ' + str(round(time.time()-t0)) + ' s') #%% check if data exist in the Pilot warehosue diff --git a/spas/acquisition.py b/spas/acquisition.py index 7701c0b..9476fbc 100644 --- a/spas/acquisition.py +++ b/spas/acquisition.py @@ -53,15 +53,18 @@ from pathlib import Path from multiprocessing import Process, Queue import shutil + import numpy as np from PIL import Image - +##### DLL for the DMD try: from ALP4 import ALP4, ALP_FIRSTFRAME, ALP_LASTFRAME from ALP4 import ALP_AVAIL_MEMORY, ALP_DEV_DYN_SYNCH_OUT1_GATE, tAlpDynSynchOutGate + # print('ALP4 is ok in Acquisition file') except: class ALP4: pass +##### DLL for the spectrometer Avantes try: from msl.equipment import EquipmentRecord, ConnectionRecord, Backend from msl.equipment.resources.avantes import MeasureCallback, Avantes @@ -69,11 +72,11 @@ class ALP4: pass from tqdm import tqdm -from .metadata import DMDParameters, MetaData, AcquisitionParameters -from .metadata import SpectrometerParameters, save_metadata, CAM, save_metadata_2arms -from .reconstruction_nn import reconstruct_process, plot_recon, ReconstructionParameters +from spas.metadata import DMDParameters, MetaData, AcquisitionParameters +from spas.metadata import SpectrometerParameters, save_metadata, CAM, save_metadata_2arms +from spas.reconstruction_nn import reconstruct_process, plot_recon, ReconstructionParameters -# Librarie for the IDS CAMERA +# DLL for the IDS CAMERA try: from pyueye import ueye, ueye_tools except: @@ -128,7 +131,7 @@ def _init_DMD() -> Tuple[ALP4, int]: # Initializing DMD dll_path = Path(__file__).parent.parent.joinpath('lib/alpV42').__str__() - + DMD = ALP4(version='4.2',libDir=dll_path) DMD.Initialize(DeviceNum=None) diff --git a/spas/metadata.py b/spas/metadata.py index 97ea932..8dfd973 100644 --- a/spas/metadata.py +++ b/spas/metadata.py @@ -21,7 +21,7 @@ import pickle ##### DLL for the DMD try: - from ALP4 import ALP4 + import ALP4 except: # in the cas the DLL of the DMD is not installed class ALP4: pass From 0abbe0025c20d57fa41669177acf862e22777aeb Mon Sep 17 00:00:00 2001 From: nducros Date: Fri, 12 May 2023 19:40:49 +0200 Subject: [PATCH 11/11] Add license and update readme file --- README.md | 52 ++++++++++------- license.md | 165 +++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 195 insertions(+), 22 deletions(-) create mode 100644 license.md diff --git a/README.md b/README.md index 000fb16..b30e96e 100644 --- a/README.md +++ b/README.md @@ -1,26 +1,18 @@ -# Single-Pixel Acquisition Software (SPAS) - Python version +# Single-Pixel Acquisition Software (SPAS) -A python toolbox for acquisition of images based on the single-pixel framework. -It has been tested using a Digital Light Processor [DLP7000](https://www.vialux.de/en/hi-speed-v-modules.html) from ViALUX GmbH and a Spectrometer [AvaSpec-ULS2048CL-EVO](https://www.avantes.com/products/spectrometers/starline/avaspec-uls2048cl-evo/) from Avantes, but may as well work for similar equipment with a few minor changes. +SPAS is python package designed for single-pixel acquisition. -## Installation (Windows only) +SPAS has been tested for controlling a [DLP7000](https://www.vialux.de/en/hi-speed-v-modules.html) Spatial Light Modulator and an [AvaSpec-ULS2048CL-EVO](https://www.avantes.com/products/spectrometers/starline/avaspec-uls2048cl-evo/) spectrometer. It should work as well for for similar equipment with a few changes. -1. Create a new environment (tested under conda) +SPAS is a companion package to the [SPyRiT](https://github.com/openspyrit/spyrit) package. -```powershell -conda create --name my_spas_env -conda activate my_spas_env -conda install -c anaconda pip -``` -2. Install the [SPyRiT](https://github.com/openspyrit/spyrit) package (tested with version 1.0.0). Typically +# Installation +The SPAS package can be installed on Linux, MacOs and Windows. However, it will be fully functional on Windows only due to DLL dependencies required for harware control. -```powershell -pip install requests torch==1.8.0+cpu torchvision==0.9.0+cpu -f https://download.pytorch.org/whl/torch_stable.html -pip install spyrit==1.0.0 -``` +We recommend using a virtual environment. -2. Clone the SPAS repository +* Clone the SPAS repository ```powershell git clone git@github.com:openspyrit/spas.git @@ -33,16 +25,15 @@ pip install -r requirements.txt pip install -e . ``` -3. Add DLLs +* Add DLLs (optional, for instrumentation control only) -The following dynamic-link libraries (DLLs) are required + The following dynamic-link libraries (DLLs) were required to control our instrumentation -* `avaspecx64.dll` provided by your Avantes distributor -* `alpV42.dll` available [here](https://www.vialux.de/en/hi-speed-download.html) by installing the entire ALP4 library + * `avaspecx64.dll` provided by your Avantes distributor + * `alpV42.dll` available [here](https://www.vialux.de/en/hi-speed-download.html) by installing the entire ALP4 library -They should be placed inside the `lib` folder -4. The typical directory structure is +* The DLLs should be placed inside the `lib` folder. The typical directory structure is ``` ├───lib @@ -62,6 +53,23 @@ They should be placed inside the `lib` folder │ ├───Cov_64x64.npy ``` +# API Documentation +https://spas.readthedocs.io/ + +# Contributors (alphabetical order) +* Thomas Baudier +* Nicolas Ducros - [Website](https://www.creatis.insa-lyon.fr/~ducros/WebPage/index.html) +* Laurent Mahieu Williame + +# How to cite? +When using SPAS in scientific publications, please cite the following paper: + +* G. Beneti-Martin, L Mahieu-Williame, T Baudier, N Ducros, "OpenSpyrit: an Ecosystem for Reproducible Single-Pixel Hyperspectral Imaging," Optics Express, Vol. 31, No. 10, (2023). https://doi.org/10.1364/OE.483937. + +# License +This project is licensed under the LGPL-3.0 license - see the [LICENSE.md](LICENSE.md) file for details + +# Getting started ## Preparation (just once) ### 1. Creating Walsh-Hadamard patterns diff --git a/license.md b/license.md new file mode 100644 index 0000000..65c5ca8 --- /dev/null +++ b/license.md @@ -0,0 +1,165 @@ + GNU LESSER GENERAL PUBLIC LICENSE + Version 3, 29 June 2007 + + Copyright (C) 2007 Free Software Foundation, Inc. + Everyone is permitted to copy and distribute verbatim copies + of this license document, but changing it is not allowed. + + + This version of the GNU Lesser General Public License incorporates +the terms and conditions of version 3 of the GNU General Public +License, supplemented by the additional permissions listed below. + + 0. Additional Definitions. + + As used herein, "this License" refers to version 3 of the GNU Lesser +General Public License, and the "GNU GPL" refers to version 3 of the GNU +General Public License. + + "The Library" refers to a covered work governed by this License, +other than an Application or a Combined Work as defined below. + + An "Application" is any work that makes use of an interface provided +by the Library, but which is not otherwise based on the Library. +Defining a subclass of a class defined by the Library is deemed a mode +of using an interface provided by the Library. + + A "Combined Work" is a work produced by combining or linking an +Application with the Library. The particular version of the Library +with which the Combined Work was made is also called the "Linked +Version". + + The "Minimal Corresponding Source" for a Combined Work means the +Corresponding Source for the Combined Work, excluding any source code +for portions of the Combined Work that, considered in isolation, are +based on the Application, and not on the Linked Version. + + The "Corresponding Application Code" for a Combined Work means the +object code and/or source code for the Application, including any data +and utility programs needed for reproducing the Combined Work from the +Application, but excluding the System Libraries of the Combined Work. + + 1. Exception to Section 3 of the GNU GPL. + + You may convey a covered work under sections 3 and 4 of this License +without being bound by section 3 of the GNU GPL. + + 2. Conveying Modified Versions. + + If you modify a copy of the Library, and, in your modifications, a +facility refers to a function or data to be supplied by an Application +that uses the facility (other than as an argument passed when the +facility is invoked), then you may convey a copy of the modified +version: + + a) under this License, provided that you make a good faith effort to + ensure that, in the event an Application does not supply the + function or data, the facility still operates, and performs + whatever part of its purpose remains meaningful, or + + b) under the GNU GPL, with none of the additional permissions of + this License applicable to that copy. + + 3. Object Code Incorporating Material from Library Header Files. + + The object code form of an Application may incorporate material from +a header file that is part of the Library. You may convey such object +code under terms of your choice, provided that, if the incorporated +material is not limited to numerical parameters, data structure +layouts and accessors, or small macros, inline functions and templates +(ten or fewer lines in length), you do both of the following: + + a) Give prominent notice with each copy of the object code that the + Library is used in it and that the Library and its use are + covered by this License. + + b) Accompany the object code with a copy of the GNU GPL and this license + document. + + 4. Combined Works. + + You may convey a Combined Work under terms of your choice that, +taken together, effectively do not restrict modification of the +portions of the Library contained in the Combined Work and reverse +engineering for debugging such modifications, if you also do each of +the following: + + a) Give prominent notice with each copy of the Combined Work that + the Library is used in it and that the Library and its use are + covered by this License. + + b) Accompany the Combined Work with a copy of the GNU GPL and this license + document. + + c) For a Combined Work that displays copyright notices during + execution, include the copyright notice for the Library among + these notices, as well as a reference directing the user to the + copies of the GNU GPL and this license document. + + d) Do one of the following: + + 0) Convey the Minimal Corresponding Source under the terms of this + License, and the Corresponding Application Code in a form + suitable for, and under terms that permit, the user to + recombine or relink the Application with a modified version of + the Linked Version to produce a modified Combined Work, in the + manner specified by section 6 of the GNU GPL for conveying + Corresponding Source. + + 1) Use a suitable shared library mechanism for linking with the + Library. A suitable mechanism is one that (a) uses at run time + a copy of the Library already present on the user's computer + system, and (b) will operate properly with a modified version + of the Library that is interface-compatible with the Linked + Version. + + e) Provide Installation Information, but only if you would otherwise + be required to provide such information under section 6 of the + GNU GPL, and only to the extent that such information is + necessary to install and execute a modified version of the + Combined Work produced by recombining or relinking the + Application with a modified version of the Linked Version. (If + you use option 4d0, the Installation Information must accompany + the Minimal Corresponding Source and Corresponding Application + Code. If you use option 4d1, you must provide the Installation + Information in the manner specified by section 6 of the GNU GPL + for conveying Corresponding Source.) + + 5. Combined Libraries. + + You may place library facilities that are a work based on the +Library side by side in a single library together with other library +facilities that are not Applications and are not covered by this +License, and convey such a combined library under terms of your +choice, if you do both of the following: + + a) Accompany the combined library with a copy of the same work based + on the Library, uncombined with any other library facilities, + conveyed under the terms of this License. + + b) Give prominent notice with the combined library that part of it + is a work based on the Library, and explaining where to find the + accompanying uncombined form of the same work. + + 6. Revised Versions of the GNU Lesser General Public License. + + The Free Software Foundation may publish revised and/or new versions +of the GNU Lesser General Public License from time to time. Such new +versions will be similar in spirit to the present version, but may +differ in detail to address new problems or concerns. + + Each version is given a distinguishing version number. If the +Library as you received it specifies that a certain numbered version +of the GNU Lesser General Public License "or any later version" +applies to it, you have the option of following the terms and +conditions either of that published version or of any later version +published by the Free Software Foundation. If the Library as you +received it does not specify a version number of the GNU Lesser +General Public License, you may choose any version of the GNU Lesser +General Public License ever published by the Free Software Foundation. + + If the Library as you received it specifies that a proxy can decide +whether future versions of the GNU Lesser General Public License shall +apply, that proxy's public statement of acceptance of any version is +permanent authorization for you to choose that version for the +Library.