Спецификация WMS позволяет серверу объединять несколько уровней в одну растровую карту. Сервер MetaCarta VMAP0 содержит множество уровней данных, таких как береговые линии, национальные границы, океан и земля. Считывание и отображение составного из нескольких слоев с сервера VMAP0. Визуализированная карта имеет пространственное разрешение 0,5 градуса на ячейку.
Найдите и обновите VMAP0 слои.
vmap0 = wmsfind('vmap0.tiles', 'SearchField', 'serverurl'); vmap0 = wmsupdate(vmap0);
Создание массива из нескольких слоев, включающих землю и океан, береговые линии и национальные границы.
layers = [refine(vmap0, 'coastline_01'); ... refine(vmap0, 'country_01'); ... refine(vmap0, 'ground_01'); ... refine(vmap0, 'inwater'); ... refine(vmap0, 'ocean')];
Извлеките составную карту. Запросите размер ячейки в 0,5 градуса, установив параметры высоты и ширины изображения. Установите значение «Transparent» равным true, чтобы все пикселы, не представляющие элементы или значения данных в слое, были установлены на прозрачное значение в полученном изображении, что позволяет создать составную карту.
[overlayImage, R] = wmsread(layers, 'Transparent', true, ... 'ImageHeight', 360, 'ImageWidth', 720);
Отображение составной карты.
figure worldmap('world') geoshow(overlayImage, R); title('Composite of VMAP0 Layers')
В этом примере показано, как можно объединить растровую карту границ с векторными данными.
Компоновка глобального растра с ячейками 1/2 градуса с помощью georefcells функция. Укажите столбцы, идущие с севера на юг, для согласования с wmsread. Результат R - объект географической растровой ссылки.
latlim = [-90 90]; lonlim = [-180 180]; cellExtent = 1/2; R = georefcells(latlim,lonlim, ... cellExtent,cellExtent,'ColumnsStartFrom','north');
Прочитайте landareas и преобразовать его в растровую карту.
land = shaperead('landareas.shp', 'UseGeoCoords', true); lat = [land.Lat]; lon = [land.Lon]; land = vec2mtx(lat, lon, zeros(R.RasterSize),R, 'filled');
Прочитайте worldrivers и преобразовать его в растровую карту.
riverLines = shaperead('worldrivers.shp','UseGeoCoords',true); rivers = vec2mtx([riverLines.Lat], [riverLines.Lon], land, R);
Объединить реки с землей.
merged = land; merged(rivers == 1) = 3;
Получение информации о системе координат из landareas shapefile с использованием shapeinfo функция. worldrivers shapefile использует ту же систему координат. Установите GeographicCRS свойства ссылочного объекта.
info = shapeinfo('landareas.shp');
R.GeographicCRS = info.CoordinateReferenceSystem;Получите изображение границ с сервера VMAP0.
vmap0 = wmsfind('vmap0.tiles', 'SearchField', 'serverurl'); vmap0 = wmsupdate(vmap0); layer = refine(vmap0, 'country_01'); height = R.RasterSize(1); width = R.RasterSize(2); [boundaries,boundariesR] = wmsread(layer, 'ImageFormat', 'image/png', ... 'ImageHeight', height, 'ImageWidth', width);
Убедитесь, что границы и объединенные растры совпадают.
isequal(boundariesR,R)
ans = logical 1
Объединить реки и земли с границами.
index = boundaries(:,:,1) ~= 255 ... & boundaries(:,:,2) ~= 255 ... & boundaries(:,:,3) ~= 255; merged(index) = 1;
Просмотрите результат.
figure worldmap(merged, R) geoshow(merged, R, 'DisplayType', 'texturemap') colormap([.45 .60 .30; 0 0 0; 0 0.5 1; 0 0 1])
Прочитайте данные о высотах и географическую ссылку на проводки для области вокруг Южного Боулдера в Колорадо. Затем обрезайте данные об отметках в меньшую область с помощью geocrop.
[fullZ,fullR] = readgeoraster('n39_w106_3arc_v2.dt1','OutputType','double'); latlim = [39.25 40.0]; lonlim = [-106 -105.5]; [Z,R] = geocrop(fullZ,fullR,latlim,lonlim);
Отображение данных фасада. Для этого создайте набор осей карты для США, постройте график данных в виде поверхности и примените соответствующую карту цветов. Просмотрите карту в 3-D, настроив положение камеры и цель. Задать вертикальное преувеличение с помощью daspectm.
figure usamap(R.LatitudeLimits, R.LongitudeLimits) geoshow(Z,R,'DisplayType','surface') demcmap(Z) title('Elevation'); cameraPosition = [218100 4367600 183700]; cameraTarget = [0 4754200 2500]; set(gca,'CameraPosition', cameraPosition, ... 'CameraTarget', cameraTarget) daspectm('m',3)
Нанесите ортоизображение на данные фасада. Для этого сначала получите названия слоев ортоизображений высокого разрешения из национальной карты USGS с помощью wmsinfo. В этом случае слой ортоизображений является единственным слоем сервера. Используйте несколько попыток подключения к серверу, если он занят.
numberOfAttempts = 5; attempt = 0; info = []; serverURL = ... 'http://basemap.nationalmap.gov/ArcGIS/services/USGSImageryOnly/MapServer/WMSServer?'; while(isempty(info)) try info = wmsinfo(serverURL); orthoLayer = info.Layer(1); catch e attempt = attempt + 1; if attempt > numberOfAttempts throw(e); else fprintf('Attempting to connect to server:\n"%s"\n', serverURL) end end end
Запрос карты слоя ортоизображений с помощью wmsread. Для отображения ортоизображений используйте geoshow и установите CData к слою.
imageHeight = size(Z,1); imageWidth = size(Z,2); orthoImage = wmsread(orthoLayer,'Latlim',R.LatitudeLimits, ... 'Lonlim',R.LongitudeLimits,'ImageHeight', imageHeight, ... 'ImageWidth', imageWidth); figure usamap(R.LatitudeLimits,R.LongitudeLimits) geoshow(Z,R,'DisplayType','surface','CData',orthoImage); title('Orthoimage Draped Over Elevation'); set(gca,'CameraPosition', cameraPosition, ... 'CameraTarget', cameraTarget) daspectm('m',3)
Файл DTED, используемый в этом примере, предоставлен Геологической службой США.