How to analyze the frequency of luminosity change of image sequence?

Hi there,

I got a set of images with different luminosity distributions, how to analyze the frequency of the luminosity change in time domain? Do we need to extract the pixel one by one to get the time series of luminosities? Besides, can we use 3D graph to show the results, with X and Y represent the location, and Z show the frequency value? Thank you very much.

I attached some images for practicing. Any help will be highly appreciated.

[Image Analyst]

Do you just want to store the mean intensity of every image into an array?

[Hao Shi]

No, I want to investigate the intensity variation at different pixels or different areas, and show their corresponding frequencies in real time with the image sequence, is it possible?

[Pippa Treacy]

did you have any success with this?

Answer 1

[Jaynik]

Hi @Hao Shi,

To analyze the frequency of luminosity changes in a set of images over time involves several steps. We need to load the sequence of images, convert them to grayscale and then calculate the difference in luminosity between the consecutive images. Following code does the same:

for i = startImage:endImage
    fileName = sprintf('%03d.tif', i); % Generates file name like '001.tif', '002.tif', etc.
    img = imread(fileName); % Read each image
    if size(img, 3) == 3
        % Image is RGB, convert to grayscale
        grayImg = rgb2gray(img);
    else
        % Image is already grayscale or in a format that does not need conversion
        grayImg = img;
    end
    
    % Assuming all images are of the same size, initialize the array if it's the first image
    if i == startImage
        [rows, cols] = size(grayImg);
        luminosityData = zeros(rows, cols, endImage - startImage + 1, 'uint8');
    end
    luminosityData(:, :, i) = grayImg; % Store grayscale image
end

% Calculate the difference in luminosity between consecutive images
luminosityChange = diff(luminosityData, 1, 3); % Calculating along the 3rd dimension (time)

% Calculate the average luminosity change
averageChange = mean(abs(luminosityChange), 3);

% Plotting the data
[X, Y] = meshgrid(1:size(averageChange, 2), 1:size(averageChange, 1));
figure;
surf(X, Y, averageChange, 'EdgeColor', 'none'); % 'EdgeColor', 'none' for a smoother appearance
xlabel('X');
ylabel('Y');
zlabel('Average Luminosity Change');
title('Luminosity Change Over Time');
colormap jet % Use a colormap that highlights changes
colorbar; % Show a color bar to interpret values

To analyze the frequency of luminosity changes, you can apply the Fourier Transform to the time series of luminosity values for each pixel. This will help in identifying the dominant frequencies of luminosity changes across the sequence of images. For each pixel, the following code also calculates the dominant frequency of luminosity change, excluding the DC component which represents the average value over time. This frequency is then stored in dominantFrequencies. A 3D surface plot is created to visualize the dominant frequencies of luminosity change across the image sequence.

% Matrix to store the dominant frequency for each pixel
dominantFrequencies = zeros(rows, cols);

% Sampling rate - assuming one image represents one time unit
Fs = 1; 

for x = 1:rows
    for y = 1:cols
        pixelSeries = squeeze(luminosityData(x, y, :)); % Extract the time series for the current pixel
        Y = fft(pixelSeries);   % Compute the Fourier Transform
        % Calculate the magnitude
        P2 = abs(Y/length(pixelSeries));
        P1 = P2(1:floor(length(pixelSeries)/2)+1);
        P1(2:end-1) = 2*P1(2:end-1);
        f = Fs*(0:(length(pixelSeries)/2))/length(pixelSeries); % Frequency domain
        
        % Finding the dominant frequency
        [~, idx] = max(P1(2:end));
        % consider frequencies > 0
        if ~isempty(idx)
            dominantFrequencies(x, y) = f(idx+1); % Adjust index for skipped DC component
        else
            dominantFrequencies(x, y) = 0;
        end
    end
end

% Visualize the data
[X, Y] = meshgrid(1:size(dominantFrequencies, 2), 1:size(dominantFrequencies, 1));
figure;
surf(X, Y, dominantFrequencies, 'EdgeColor', 'none');
xlabel('X');
ylabel('Y');
zlabel('Dominant Frequency');
title('Dominant Frequency of Luminosity Change Over Time');
colormap jet % Use a colormap that highlights changes
colorbar; % Show a color bar to interpret values

To investigate the intensity variation at different pixels or areas and display their corresponding frequencies in real time with the image sequence, you would need a slightly different approach than what was previously discussed. The earlier code calculates the dominant frequency of luminosity change for each pixel across the entire image sequence and visualizes it as a static 3D plot. To analyze and visualize the intensity variation and corresponding frequencies in real time or dynamically with the image sequence, you would need to incorporate interactive or animated visualizations. The following code create an animation that iterates through each time frame, showing the original image and its frequency analysis.

figure;
for t = 1:size(luminosityData, 3)
    subplot(1, 2, 1);
    imshow(luminosityData(:, :, t));
    title(sprintf('Original Image at Time %d', t));
    
    subplot(1, 2, 2);
    imagesc(dominantFrequencies);
    title('Frequency Analysis');
    colormap jet;
    colorbar;
    
    pause(0.5); % Adjust the pause for desired animation speed
end

You can read the following documentations to read more about these functions:

• meshgrid: https://www.mathworks.com/help/matlab/ref/meshgrid.html
• squeeze: https://www.mathworks.com/help/matlab/ref/squeeze.html
• fft: https://https.mathworks.com/help/matlab/ref/fft.html
• imagesc: https://www.mathworks.com/help/matlab/ref/imagesc.html

Hope this helps!