I want to plot in a different barcolors

Dear All,

I have this code and when I run it it shows both set of data in jet colorbar while I want to distinguish these two from eachother? Would you please advice me what to do?

Best,

Ara

% Initialize video writer for combined visualization
video_filename = 'ROTI_S4_Map_Visualization_MovingCircles_4.avi';
video_writer = VideoWriter(video_filename, 'Motion JPEG AVI');
video_writer.FrameRate = 30;
open(video_writer);

% Define the date range for filtering
start_date = datetime(2023, 4, 23, 0, 0, 0);
end_date = datetime(2023, 4, 23, 23, 59, 59);

% Initialize figure and map
figure('Position', [100, 100, 1200, 800]);
latlim = [35 70];
lonlim = [-10 40];

% Pre-read S4 data for all stations
station_data = struct();
for station_num = 1:29
    filename_s4 = sprintf('Station%d.xlsx', station_num);
    try
        data_s4 = readmatrix(filename_s4, 'OutputType', 'string');
        station_data(station_num).data = data_s4;
    catch
        disp(['File ', filename_s4, ' not found. Skipping.']);
    end
end

% Loop over each 5-minute interval
fileDates = {'20230423'};
for d = 1:length(fileDates)
    date_str = fileDates{d};
    time_bins = start_date:minutes(5):end_date;

    for t = 1:length(time_bins)-1
        % Clear current frame
        cla;

        % Re-draw the base map
        ax = worldmap(latlim, lonlim);
        setm(ax, 'FFaceColor', [0.5 0.7 0.9]);
        load coastlines;
        geoshow(coastlat, coastlon, 'DisplayType', 'polygon', 'FaceColor', [0.6 0.8 0.6]);
        geoshow('landareas.shp', 'FaceColor', [0.8 0.8 0.8]);
        hold on;

        % Define current time window
        current_time = time_bins(t);
        hour_str = datestr(current_time, 'HH');
        minute_str = datestr(current_time, 'MM');
        filename_vtec = sprintf('C:\\Users\\I7-2600K\\Desktop\\data\\%s_%s%s00_roti.dat', ...
                                 date_str, hour_str, minute_str);

        if isfile(filename_vtec)
            data_ROTI = readmatrix(filename_vtec);
            if size(data_ROTI, 2) >= 3
                lon_ROTI = data_ROTI(:, 1);
                lat_ROTI = data_ROTI(:, 2);
                ROTI = data_ROTI(:, 3);
                scatterm(lat_ROTI, lon_ROTI, 30, ROTI, 'filled');
                colormap(ax, jet);
                colorbar;
                caxis([0 0.5]);
                cb_roti = colorbar('southoutside');
                cb_roti.Label.String = 'ROTI';
                title(['ROTI Map - ', datestr(current_time)]);
            end
        else
            disp(['ROTI file not found: ', filename_vtec]);
        end

        % Plot S4 as moving circles
        for station_num = 1:29
            if isfield(station_data(station_num), 'data')
                data_s4 = station_data(station_num).data;
                time_str_s4 = data_s4(:, 1);
                S4 = str2double(data_s4(:, 3));
                Lon_IPP = str2double(data_s4(:, 15));
                Lat_IPP = str2double(data_s4(:, 16));
                elevation = str2double(data_s4(:, 5));
                time_dt_s4 = datetime(time_str_s4, 'InputFormat', 'yyyy-MM-dd''T''HH:mm:ss.SSS');

                % Filter for current time interval and elevation > 30 deg
                valid_idx = (time_dt_s4 >= time_bins(t)) & ...
                            (time_dt_s4 < time_bins(t+1)) & ...
                            (elevation >= 30);

                if any(valid_idx)
                    S4_filtered = S4(valid_idx);
                    Lon_IPP_filtered = Lon_IPP(valid_idx);
                    Lat_IPP_filtered = Lat_IPP(valid_idx);

                    for i = 1:length(S4_filtered)
                        color = [0 1 0]; % Green for low scintillation
                        if S4_filtered(i) > 0.05
                            color = [1 0 0]; % Red for high scintillation
                        end
                        scatterm(Lat_IPP_filtered(i), Lon_IPP_filtered(i), ...
                                 100, color, 'o', 'filled', ...
                                 'MarkerEdgeColor', 'k', 'LineWidth', 0.5);
      
                    end
                end
            end
        end

        colormap(ax, parula); % Can be changed
        caxis([0 0.3]);
        cb_s4 = colorbar('eastoutside');
        cb_s4.Label.String = 'S4';


        drawnow;
        frame = getframe(gcf);
        writeVideo(video_writer, frame);
    end
end

close(video_writer);
disp(['Video saved to ', video_filename]);

Answer 1

[Sourabh]

Hey @Ara

The issue you are facing arises because you are using two datasets (ROTI and S4) but applying only one colormap at a time and overwriting it during plotting. This causes confusion in the colorbar because the last “colormap()” and “colorbar()” command wins, and both datasets might appear as though they share the same scale and colors, even though they represent different phenomena.

You need to create two separate colorbars, each corresponding only to its dataset. Right now, both datasets are plotted on the same axes (ax), so they share a single colormap context. To separate them visually kindly refer to the following steps:

• Keep your main plotting on a single axes (ax).
• Use invisible dummy axes to separately handle each colormap + colorbar:

ax_roti for ROTI using jet

ax_s4 for S4 using parula

• Plot ROTI and S4 as usual, but apply colormaps only to their respective dummy axes.

Modified Code Snippet:

Replace your ROTI + S4 plotting and colorbar code with this:

% ----------- ROTI Plot ----------- 
scatterm(lat_ROTI, lon_ROTI, 30, ROTI, 'filled'); 


% Create dummy axes for ROTI colormap and colorbar 
ax_roti = axes('Position', ax.Position, 'Color', 'none', 'Visible', 'off'); 
colormap(ax_roti, jet); 
cb_roti = colorbar(ax_roti, 'southoutside'); 
cb_roti.Label.String = 'ROTI'; 
cb_roti.Position = [0.35 0.08 0.3 0.02]; % Adjust as needed 
caxis(ax_roti, [0 0.5]); 


% ----------- S4 Plot ----------- 
for station_num = 1:29 
    if isfield(station_data(station_num), 'data') 
        … 
        … 
        … 
        end 
    end 
end 
 

% Create dummy axes for S4 colormap and colorbar 
ax_s4 = axes('Position', ax.Position, 'Color', 'none', 'Visible', 'off'); 
colormap(ax_s4, parula); 
cb_s4 = colorbar(ax_s4, 'eastoutside'); 
cb_s4.Label.String = 'S4'; 
cb_s4.Position = [0.88 0.25 0.02 0.5]; % Adjust as needed to avoid overlap
caxis(ax_s4, [0 0.3]); 

For more information on “colormap” and “colorbar”, kindly refer the following MATLAB documentations:

• https://www.mathworks.com/help/matlab/ref/colormap.html
• https://www.mathworks.com/help/matlab/ref/colorbar.html

[Ara]

Dear Sourabh,

Thank you very much for responsing my question, and corrected my code. I replaced it as you suggested but the Matlab getting busy without shows the plots. Do you have any idea of what to do?

Best regards,

Ara

[Image Analyst]

@Ara can you post a screenshot?

If you have any more questions, then attach your data and code to read it in with the paperclip icon after you read this:

TUTORIAL: How to ask a question (on Answers) and get a fast answer

[Ara]

Dear Image Analysis,

Here is the screenshot. Please see attached.

Unfortunately, I am not able to share the data.