-
Notifications
You must be signed in to change notification settings - Fork 36
/
CatPoly2GRIDobj.m
105 lines (93 loc) · 3.43 KB
/
CatPoly2GRIDobj.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
function [OUT,look_table]=CatPoly2GRIDobj(DEM,poly_shape,field,varargin)
%
% Usage:
% [GRIDobj,look_up_table]=CatPoly2GRIDobj(DEM,poly_shape,field);
%
% Description:
% Function to convert a categorical polygon shape file (e.g. a digitzed geologic map) to a GRIDobj. Can
% be useful for use in 'ProcessRiverBasins'
%
% Required Inputs:
% DEM - DEM that you want the output to match
% poly_shape - name or path to shapefile containing the categorical data
% field - field name of categorical data within the shapefile
%
% Optional Input:
% table_in - Optional input to manually provide a table to use as the look_table. This table must contain
% two columns where the first column is n x 1 array of scalar integers and the second column is a n x 1
% cell array containing the category names. These categories must match entries in the 'field' of the
% provided 'poly_shape' or all of the values in the 'OUT' will be 0.
%
%
% Outputs:
% OUT - GRIDobj of the same size as DEM where values correspond to categorical data
% as defined in the look_table
% look_table - nx2 table with columns Numbers and Categories that serves as a lookup table
% to convert between the numbers and the original categories.
%
% Example:
% [GEO,geo_table]=CatPoly2GRIDobj(DEM,'geologic_map.shp','rtype');
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Function Written by Adam M. Forte - Updated : 07/05/19 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Parse Inputs
p = inputParser;
p.FunctionName = 'CatPoly2GRIDobj';
addRequired(p,'DEM',@(x) isa(x,'GRIDobj'));
addRequired(p,'poly_shape',@(x) ischar(x));
addRequired(p,'field',@(x) ischar(x));
addParameter(p,'table_in',[],@(x) istable(x));
parse(p,DEM,poly_shape,field,varargin{:});
DEM=p.Results.DEM;
poly_shape=p.Results.poly_shape;
field=p.Results.field;
table_in=p.Results.table_in;
% Validate table if provided
if ~isempty(table_in)
if size(table_in,2)~=2
error('Provided "table_in" must have two columns');
elseif table_in.(1)(1)~=0
error('First entry of first column of "table_in" must equal 0');
elseif ~iscell(table_in.(2))
error('Second column of "table_in" must be a cell array');
elseif ~ischar(table_in.(2){1})
error('Entries in second column of "table_in" must be characters');
end
end
% Read in shape and covert to a table
PS=shaperead(poly_shape);
TS=struct2table(PS);
% Separate out the field of interest
Foi=TS.(field);
if isempty(table_in)
% Generate unique list and the output lookup table
Categories=unique(Foi);
Numbers=[1:numel(Categories)]; Numbers=Numbers';
% Add an 'undef' category to deal with zeros that appear because of read errors
Numbers=vertcat(0,Numbers);
Categories=vertcat('undef',Categories);
look_table=table(Numbers,Categories);
else
% Load provided table into specific table format
Categories=table_in.(2);
Numbers=table_in.(1);
look_table=table(Numbers,Categories);
end
% Replace categorical with number
for ii=1:numel(PS)
Eoi=PS(ii,1).(field);
ix=find(strcmp(Categories,Eoi));
PS(ii,1).replace_number=Numbers(ix);
end
% Run polygon2GRIDobj with pad to avoid ring of 0s bug
DEMp=GRIDobj(DEM,'logical');
DEMp.Z(:,:)=true;
DEMp=pad(DEMp,1,false);
[OUT]=polygon2GRIDobj(DEMp,PS,'replace_number');
OUT=crop(OUT,DEMp);
% Remove nonexistent category-number pairs from look_table
pres=unique(OUT.Z(:));
ix=ismember(look_table.Numbers,pres);
look_table=look_table(ix,:);
end