Merge Rasters of Cells

This example shows how to import two adjacent raster data files made of cells, merge the data into a single raster, and display the merged raster on a map.

The files used in this example, concord_ortho_e.tif and concord_ortho_w.tif, contain east-west adjacent images with concord_ortho_e.tif to the east of concord_ortho_w.tif. The cell extents and world y-limits of the images are identical.

Import the images using the readgeoraster function. The outputs eastA and westA are arrays that contain the image data, and the outputs westR and eastR are MapCellsReference objects that contain spatial referencing information.

[eastA,eastR] = readgeoraster('concord_ortho_e.tif');
[westA,westR] = readgeoraster('concord_ortho_w.tif');

Create a merged raster by combining the arrays. The images are adjacent and are made up of cells, so the eastern boundary of westA aligns with the western boundary of eastA.

mergedA = [westA eastA];

Spatially reference the merged raster by creating a raster reference object. You can create a reference object for a projected raster of cells by using the maprefcells function.

Specify the world x- and y-limits of the raster. The x-limits of the merged raster are the minimum x-limit of the western raster and the maximum x-limit of the eastern raster. The world y-limits of the merged raster are the same as the y-limits of the imported rasters. Create the reference object.

xlimits = [westR.XWorldLimits(1) eastR.XWorldLimits(2)];
ylimits = westR.YWorldLimits;
mergedR = maprefcells(xlimits,ylimits,size(mergedA));

The columns of arrays imported using the readgeoraster function start from the north. Therefore, set the ColumnsStartFrom property of the reference object to 'north'.

mergedR.ColumnsStartFrom = 'north';

Display the merged image on a map by using the mapshow function.

mapshow(mergedA,mergedR,'DisplayType','image')

Merge Rasters of Postings

This example shows how to import two adjacent raster data files made of posting points, merge the data into a single raster, and display the merged raster on a map.

The files used in this example, n39_w106_3arc_v2.dt1 and n40_w106_3arc_v2.dt1, contain north-south adjacent elevation grids with n39_w106_3arc_v2.dt1 to the south of n40_w106_3arc_v2.dt1. The northern latitude limit of the grid in n39_w106_3arc_v2.dt1 is the same as the southern latitude limit of the grid in n40_w106_3arc_v2.dt1. The spacing of posting points and the longitude limits of the grids are identical.

Import the grids using the readgeoraster function. The outputs southZ and northZ are arrays of type double that contain the elevation data, and the outputs southR and northR are GeographicPostingsReference objects that contain spatial referencing information.

[southZ,southR] = readgeoraster('n39_w106_3arc_v2.dt1','OutputType','double');
[northZ,northR] = readgeoraster('n40_w106_3arc_v2.dt1','OutputType','double');

Create a merged raster by combining the arrays. The boundaries of posting point rasters are made up of the outermost posting points, and the columns of arrays imported using the readgeoraster function start from the north. Therefore, the southernmost row of northZ and the northernmost row of southZ are the same. To avoid a duplicate row in the merged raster, remove the southernmost row of northZ before combining the arrays.

northZ(end,:) = [];
mergedZ = [northZ; southZ];

Spatially reference the merged raster by creating a raster reference object. You can create a reference object for a geographic raster of posting points by using the georefpostings function.

Specify the latitude and longitude limits of the raster. The latitude limits are the minimum latitude of the southern raster and the maximum latitude of the northern raster. The longitude limits of the merged raster are the same as the longitude limits of the imported rasters. Create the reference object.

latlim = [southR.LatitudeLimits(1) northR.LatitudeLimits(2)];
lonlim = southR.LongitudeLimits;
mergedR = georefpostings(latlim,lonlim,size(mergedZ));

Set the ColumnsStartFrom and GeographicCRS properties of the reference object so that they match the properties of the imported reference objects.

mergedR.ColumnsStartFrom = southR.ColumnsStartFrom;
mergedR.GeographicCRS = southR.GeographicCRS;

Display the merged raster on a map. Create a axesm-based map by specifying the latitude and longitude limits of the data. Then, display the data as a surface by using the geoshow function. Apply a colormap appropriate for elevation data by using the demcmap function.

usamap(mergedR.LatitudeLimits,mergedR.LongitudeLimits)
geoshow(mergedZ,mergedR,'DisplayType','surface')
demcmap(mergedZ)

The elevation data used in this example is from the U.S. Geological Survey.