Calculations/Data-Analyzer/+Scripts/BECToStripesToDroplets/runAndCompareAngularSpectralDistributionAnalysis.m

191 lines
7.7 KiB
Matlab
Raw Permalink Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

%% ===== Stripes-Droplets: Compare Mean and Variance Across Data Sources =====
% Author: Karthik
% Date: 2025-10-29
% Version: 1.2
%
% Description:
% Loops through multiple data sources, performs angular spectral analysis,
% extracts the first two cumulants (Mean, Variance), and plots them
% together in a 1×2 tiled layout
% Each curve can be labeled manually.
%
% Saves both figure (.fig) and data (.mat) for replotting later.
clear; clc; close all;
%% ===== Define Data Sources and Labels =====
dataSourcesList = {
struct('sequence', 'StructuralPhaseTransition', 'date', '2025/10/13', 'runs', [11]), ...
struct('sequence', 'StructuralPhaseTransition', 'date', '2025/10/13', 'runs', [12]), ...
struct('sequence', 'StructuralPhaseTransition', 'date', '2025/10/14', 'runs', [0]), ...
struct('sequence', 'StructuralPhaseTransition', 'date', '2025/10/14', 'runs', [1]), ...
struct('sequence', 'StructuralPhaseTransition', 'date', '2025/10/14', 'runs', [2]), ...
struct('sequence', 'StructuralPhaseTransition', 'date', '2025/10/14', 'runs', [3]), ...
struct('sequence', 'StructuralPhaseTransition', 'date', '2025/10/14', 'runs', [4])
};
curveLabels = {'B = 2.55 G', 'B = 2.50 G', 'B = 2.45 G', ...
'B = 2.40 G', 'B = 2.35 G', 'B = 2.30 G', ...
'B = 2.20 G'};
%% ===== Analysis Options =====
options = struct();
% File paths
options.baseDataFolder = '//DyLabNAS/Data';
options.FullODImagesFolder = 'E:/Data - Experiment/FullODImages/202510';
options.measurementName = 'StripesToDroplets';
scriptFullPath = mfilename('fullpath');
options.saveDirectory = fileparts(scriptFullPath);
% Camera / imaging settings
options.cam = 4; % 1 - ODT_1_Axis_Camera; 2 - ODT_2_Axis_Camera; 3 - Horizontal_Axis_Camera;, 4 - Vertical_Axis_Camera;
options.angle = 0; % angle by which image will be rotated
options.center = [1420, 2040];
options.span = [200, 200];
options.fraction = [0.1, 0.1];
options.pixel_size = 5.86e-6; % in meters
options.magnification = 23.2;
options.ImagingMode = 'HighIntensity';
options.PulseDuration = 5e-6; % in s
% Fourier analysis settings
options.theta_min = deg2rad(0);
options.theta_max = deg2rad(180);
options.N_radial_bins = 500;
options.Radial_Sigma = 2;
options.Radial_WindowSize = 5; % odd number
options.k_min = 1.0500; % μm⁻¹
options.k_max = 2.7000; % μm⁻¹
options.N_angular_bins = 360;
options.Angular_Threshold = 75;
options.Angular_Sigma = 2;
options.Angular_WindowSize = 5;
options.zoom_size = 50;
% Processing flags
options.skipUnshuffling = false;
options.skipNormalization = false;
options.skipFringeRemoval = true;
options.skipPreprocessing = true;
options.skipMasking = true;
options.skipIntensityThresholding = true;
options.skipBinarization = true;
options.skipFullODImagesFolderUse = false;
options.skipSaveData = false;
options.skipSaveFigures = true;
options.skipSaveProcessedOD = true;
options.skipLivePlot = false;
options.showProgressBar = true;
options.font = 'Bahnschrift';
% ===== Initialize container =====
cumulantData = struct(); % will be indexed as cumulantData(thetaIdx, dataSourceIdx)
% ===== Specify theta values to analyze =====
desiredTheta = [pi/6, pi/3]; % you can add more angles here
%% ===== Loop over data sources =====
for d = 1:numel(dataSourcesList)
dataSources = {dataSourcesList{d}};
[options.selectedPath, options.folderPath] = Helper.selectDataSourcePath(dataSources, options);
% --- Read all images for this data source once ---
[od_imgs, scan_parameter_values, scan_reference_values, file_list] = ...
Helper.collectODImages(options);
spectral_analysis_results = Analyzer.extractFullAngularSpectralDistribution(od_imgs, options);
results = Analyzer.recenterSpectralCurves( ...
spectral_analysis_results.S_theta_norm_all, ...
spectral_analysis_results.theta_vals/pi, ...
scan_reference_values, ...
'SearchRange', [0 90]);
mask = results.x_values >= 0 & results.x_values <= 1;
results.x_values = results.x_values(mask);
for i = 1:numel(results.curves)
results.curves{i} = results.curves{i}(:, mask);
end
% --- Loop over theta values for this data source ---
for t = 1:numel(desiredTheta)
th = desiredTheta(t);
thetaVals = results.x_values * pi;
[~, thIdx] = min(abs(thetaVals - th));
N_params = numel(results.curves);
meanVals = zeros(1, N_params);
varVals = zeros(1, N_params);
for i = 1:N_params
reps_values = results.curves{i}(:, thIdx);
kappa = Calculator.computeCumulants(reps_values, 2);
meanVals(i) = kappa(1);
varVals(i) = kappa(2);
end
% --- Store cumulants ---
cumulantData(t,d).theta = th;
cumulantData(t,d).label = curveLabels{d};
cumulantData(t,d).mean = meanVals;
cumulantData(t,d).variance = varVals;
cumulantData(t,d).scan_vals = scan_reference_values;
end
end
%% ===== Plot Mean and Variance per theta =====
N_theta = size(cumulantData,1); % number of theta values
N_data = size(cumulantData,2); % number of data sources
% --- Colormap (trimmed coolwarm) ---
fullCmap = Colormaps.coolwarm(256);
blueSide = fullCmap(1:100,:);
redSide = fullCmap(157:end,:);
trimmedCmap = [blueSide; redSide];
cmap = trimmedCmap(round(linspace(1,size(trimmedCmap,1),N_data)), :);
for t = 1:N_theta
fig = figure('Color','w','Position',[100 100 1200 550]);
tName = sprintf('\\theta = %.0f°', round(rad2deg(cumulantData(t,1).theta)));
sgtitle(['Cumulants at ', tName], 'FontName', options.font, 'FontSize', 18);
tLayout = tiledlayout(1,2,'TileSpacing','Compact','Padding','Compact');
% --- Mean (κ₁) ---
ax1 = nexttile; hold(ax1,'on');
for d = 1:N_data
plot(ax1, cumulantData(t,d).scan_vals, cumulantData(t,d).mean, '-o', ...
'Color', cmap(d,:), 'LineWidth', 2, 'MarkerSize', 8, ...
'MarkerFaceColor', cmap(d,:), 'DisplayName', cumulantData(t,d).label);
end
xlabel(ax1, '\alpha (degrees)', 'FontName', options.font, 'FontSize', 14);
ylabel(ax1, '\kappa_1', 'FontName', options.font, 'FontSize', 14);
title(ax1, 'Mean', 'FontName', options.font, 'FontSize', 16);
grid(ax1,'on'); legend(ax1,'Location','northwest','FontSize',12);
set(ax1,'FontName',options.font,'FontSize',14);
% --- Variance (κ₂) ---
ax2 = nexttile; hold(ax2,'on');
for d = 1:N_data
plot(ax2, cumulantData(t,d).scan_vals, cumulantData(t,d).variance, '-o', ...
'Color', cmap(d,:), 'LineWidth', 2, 'MarkerSize', 8, ...
'MarkerFaceColor', cmap(d,:), 'DisplayName', cumulantData(t,d).label);
end
xlabel(ax2, '\alpha (degrees)', 'FontName', options.font, 'FontSize', 14);
ylabel(ax2, '\kappa_2', 'FontName', options.font, 'FontSize', 14);
title(ax2, 'Variance', 'FontName', options.font, 'FontSize', 16);
grid(ax2,'on'); legend(ax2,'Location','northwest','FontSize',12);
set(ax2,'FontName',options.font,'FontSize',14);
% --- Save figure ---
saveFileNameFig = fullfile(options.saveDirectory, ...
sprintf('StripesToDroplets_MeanVariance_theta_%d.fig', round(rad2deg(cumulantData(t,1).theta))));
savefig(fig, saveFileNameFig);
end
% --- Save all cumulant data ---
saveFileNameMat = fullfile(options.saveDirectory, 'StripesToDroplets_MeanVariance_AllThetaData.mat');
save(saveFileNameMat, 'cumulantData', 'options');