Approach for plotting with an indeterminate number of colormaps

2 次查看(过去 30 天)
The following code allows multiple landsat tiles to be plotted using the m_map package. The variable IMAGERY_PATH is a cell which is intended to hold an indeterminate number of paths in which the landsat band .tif files reside. Presently, only the first (or last?) tile is rendered with a sensible color. It would seem that a new colormap should be used for each tile.
Any comments on suggested approach are greatly appreciated. I'm happy to provide a link to a set of relatively small data files.
% plot_multiple_landsat_tiles.m
addpath 'C:\matlab\L8read\'
addpath 'C:\matlab\m_map1_4\m_map'
% constants
proj = 'UTM';
datum = 'wsg84';
%%{
LON = [-81.5, -74.5];
LAT = [-13.5, -8.5];
%}
%{
LON = [-80, -77];
LAT = [-13.5, -12.5];
%}
% landsat data
% IMAGERY_PATH = {'...\LC08_L1TP_007069_20170408_20170414_01_T1\'; ...
% '...\LC08_L1TP_009067_20170406_20170414_01_T1\'};
% IMAGERY_PATH = {'...\LC08_L1TP_007069_20170408_20170414_01_T1\'; ...
% '...\LC08_L1TP_006071_20170401_20170414_01_T1\'};
IMAGERY_PATH = {'...\LC08_L1TP_007069_20170408_20170414_01_T1\'};
num_source_images = length(IMAGERY_PATH);
% instantiate figure before loop
ax1 = axes; % instantiate axes
hold(ax1,'on') % hold on
set(gcf,'Color','w')
set(ax1,'Color','none')
set(ax1,'Box','on') % box on
m_proj(proj, ...
'long', LON, ...
'lat', LAT, ...
'ellipse', datum); % map projection
m_grid; % show grid
set(ax1, 'Color', 'none') % set axis background color to none
m_gshhs_f('patch', [0.8 0.8 0.8]); % coastline patch
for i = 1: num_source_images
imagery_path = IMAGERY_PATH{i};
[I_B2, ~, ~, info] = L8read(imagery_path, LAT, LON, 'band', 2, 'stretch');
[I_B3, ~, ~, ~] = L8read(imagery_path, LAT, LON, 'band', 3, 'stretch');
[I_B4, ~, ~, ~] = L8read(imagery_path, LAT, LON, 'band', 4, 'stretch');
% composite
I = cat(3, I_B4, I_B3, I_B2);
% convert the image to an indexed image
n = 65536;
[I_indexed, color_map] = rgb2ind(I, n);
% for support of nan's, convert to type double
I_double = im2double(I_indexed, 'indexed');
% where the image pixel is 1, let it be nan
I_double(I_double == 1) = NaN;
% get the corner longitude, lattitude coordinates
image_corner_lon = [info.CornerCoords.Lon(1,1) info.CornerCoords.Lon(1,2)];
image_corner_lat = [info.CornerCoords.Lat(1,1) info.CornerCoords.Lat(1,3)];
% image corners
[image_corner_mmapx, image_corner_mmapy] = m_ll2xy(image_corner_lon, image_corner_lat);
% define colormap
colormap(ax1, color_map);
% display the composite landsat imagery
im = image(image_corner_mmapx, image_corner_mmapy, I_double);
% not sure what this does, if anything
im.AlphaData = double(~isnan(I_double));
clear I_B4 I_B3 I_B2 x y info I I_double
end

回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Convert Image Type 的更多信息

产品


版本

R2018b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by