From d4c0dc3f07f4096de088380cedd0ece126fd6263 Mon Sep 17 00:00:00 2001 From: Karthik Chandrashekara Date: Fri, 11 Apr 2025 03:02:00 +0200 Subject: [PATCH] Integrated S-theta, k-means clustering code added. --- Data-Analyzer/fourierAnalysis.m | 109 ++++++++++++++++++++++++++++---- 1 file changed, 98 insertions(+), 11 deletions(-) diff --git a/Data-Analyzer/fourierAnalysis.m b/Data-Analyzer/fourierAnalysis.m index 3f05c5c..a29b3d3 100644 --- a/Data-Analyzer/fourierAnalysis.m +++ b/Data-Analyzer/fourierAnalysis.m @@ -80,6 +80,7 @@ end %% Run Fourier analysis over images fft_imgs = cell(1, nimgs); +integrated_sf = zeros(1, nimgs); % Create VideoWriter object for movie videoFile = VideoWriter('Single_Shot_FFT.mp4', 'MPEG-4'); @@ -92,7 +93,7 @@ for k = 1 : length(od_imgs) IMG = od_imgs{k}; [IMGFFT, IMGBIN] = computeFourierTransform(IMG); - figure(2); + figure(1); clf set(gcf,'Position',[500 100 1000 800]) t = tiledlayout(2, 2, 'TileSpacing', 'compact', 'Padding', 'compact'); % 1x4 grid @@ -146,8 +147,8 @@ for k = 1 : length(od_imgs) hYLabel = ylabel('Y (pixels)', 'Interpreter', 'tex'); hTitle = title('Denoised - Masked - Binarized', 'Interpreter', 'tex'); - set([hXLabel, hYLabel, hL], 'FontName', font) - set([hXLabel, hYLabel, hL], 'FontSize', 14) + set([hXLabel, hYLabel], 'FontName', font) + set([hXLabel, hYLabel], 'FontSize', 14) set(hTitle, 'FontName', font, 'FontSize', 16, 'FontWeight', 'bold'); % Set font and size for title ax3 = nexttile; @@ -172,11 +173,11 @@ for k = 1 : length(od_imgs) hXLabel = xlabel('X (pixels)', 'Interpreter', 'tex'); hYLabel = ylabel('Y (pixels)', 'Interpreter', 'tex'); hTitle = title('Fourier Power Spectrum', 'Interpreter', 'tex'); - set([hXLabel, hYLabel, hL, hText], 'FontName', font) - set([hXLabel, hYLabel, hL], 'FontSize', 14) + set([hXLabel, hYLabel, hText], 'FontName', font) + set([hXLabel, hYLabel], 'FontSize', 14) set(hTitle, 'FontName', font, 'FontSize', 16, 'FontWeight', 'bold'); % Set font and size for title - % Plot the angular structure factor + % Plot the angular distribution %{ nexttile [theta_vals, angular_intensity] = computeAngularDistribution(zoomedIMGFFT, 10, 20, 100, 75); @@ -190,19 +191,17 @@ for k = 1 : length(od_imgs) % Plot the angular structure factor nexttile [theta_vals, S_theta] = computeAngularStructureFactor(zoomedIMGFFT, 10, 20, 180, 75, 2); - + integrated_sf(k) = trapz(theta_vals, S_theta); plot(theta_vals/pi, S_theta,'Linewidth',2); set(gca, 'FontSize', 14); % For tick labels only hXLabel = xlabel('\theta (\pi)', 'Interpreter', 'tex'); hYLabel = ylabel('S(\theta)', 'Interpreter', 'tex'); hTitle = title('Angular Structure Factor', 'Interpreter', 'tex'); - set([hXLabel, hYLabel, hL, hText], 'FontName', font) - set([hXLabel, hYLabel, hL], 'FontSize', 14) + set([hXLabel, hYLabel, hText], 'FontName', font) + set([hXLabel, hYLabel], 'FontSize', 14) set(hTitle, 'FontName', font, 'FontSize', 16, 'FontWeight', 'bold'); % Set font and size for title grid on - - drawnow pause(0.5) @@ -214,6 +213,94 @@ end % Close the video file close(videoFile); +%% Track integrated structure factor across the transition + +% Assuming theta_values and integrated_sf are column vectors (or row vectors of same length) +[unique_theta, ~, idx] = unique(theta_values); + +% Preallocate arrays +mean_sf = zeros(size(unique_theta)); +stderr_sf = zeros(size(unique_theta)); + +% Loop through each unique theta and compute mean and standard error +for i = 1:length(unique_theta) + group_vals = integrated_sf(idx == i); + mean_sf(i) = mean(group_vals); + stderr_sf(i) = std(group_vals) / sqrt(length(group_vals)); % standard error = std / sqrt(N) +end + +figure(2); +set(gcf,'Position',[100 100 950 750]) +errorbar(unique_theta, mean_sf, stderr_sf, 'o--', ... + 'LineWidth', 1.5, 'MarkerSize', 6, 'CapSize', 5); +set(gca, 'FontSize', 14); % For tick labels only +hXLabel = xlabel('\alpha (degrees)', 'Interpreter', 'tex'); +hYLabel = ylabel('Integrated S(\theta)', 'Interpreter', 'tex'); +hTitle = title('Change during rotation', 'Interpreter', 'tex'); +set([hXLabel, hYLabel], 'FontName', font) +set([hXLabel, hYLabel], 'FontSize', 14) +set(hTitle, 'FontName', font, 'FontSize', 16, 'FontWeight', 'bold'); % Set font and size for title +grid on + +%% k-means Clustering + +% Reshape to column vector +X = mean_sf(:); + +% Determine the number of clusters to try (you can experiment with different values here) +optimalClusters = 3; % Based on prior knowledge or experimentation + +% Set the random seed for reproducibility +rng(42); + +% Specify initialization method ('plus' is the default) +startMethod = 'plus'; % Options: 'uniform', 'plus', 'sample' + +% Apply K-means clustering with controlled initialization +[idx, C] = kmeans(X, optimalClusters, 'Start', startMethod); + +% Plot the results +figure(3); +set(gcf,'Position',[100 100 950 750]) +hold on; + +% Plot error bars with mean_sf and stderr_sf +errorbar(unique_theta, mean_sf, stderr_sf, 'o--', ... + 'LineWidth', 1.5, 'MarkerSize', 6, 'CapSize', 5); + +% Scatter plot for data points (showing clusters) +scatter(unique_theta, X, 50, idx, 'filled'); + +% Get the current y-axis limits +current_ylim = ylim; + +% Generate colors for each cluster +colors = lines(optimalClusters); % Create distinct colors for each cluster + +% Loop through each cluster and fill the regions +for i = 1:optimalClusters + % Find indices of data points that belong to the current cluster + clusterIdx = find(idx == i); + + % Find the range of x-values for this cluster + x_min = unique_theta(clusterIdx(1)); % Starting x-value for the cluster + x_max = unique_theta(clusterIdx(end)); % Ending x-value for the cluster + + % Fill the region corresponding to the cluster + fill([x_min, x_max, x_max, x_min], ... + [current_ylim(1), current_ylim(1), current_ylim(2), current_ylim(2)], ... + colors(i, :), 'EdgeColor', 'none', 'FaceAlpha', 0.3); % Add transparency +end + +hXLabel = xlabel('\alpha (degrees)', 'Interpreter', 'tex'); +hYLabel = ylabel('Integrated S(\theta)', 'Interpreter', 'tex'); +hTitle = title('Regimes identified by k-means clustering', 'Interpreter', 'tex'); +set([hXLabel, hYLabel], 'FontName', font) +set([hXLabel, hYLabel], 'FontSize', 14) +set(hTitle, 'FontName', font, 'FontSize', 16, 'FontWeight', 'bold'); % Set font and size for title +grid on; +hold off; + %% Helper Functions function [IMGFFT, IMGBIN] = computeFourierTransform(I) % computeFourierSpectrum - Computes the 2D Fourier power spectrum