forked from gustavo-gomes-ghg/My_Matlab
-
Notifications
You must be signed in to change notification settings - Fork 2
/
xyz2grd.m
56 lines (49 loc) · 1.56 KB
/
xyz2grd.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
function [lat,lon,z] = xyz2grd(xyz,fname)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% function [lat,lon,z] = xyz2grd(xyz,fname)
%
% Takes xyz geographic file and returns a lat, lon and z array
% Optionally creates a .grd and .mask file from the z array.
%
% INPUT
% xyz - (2D Array) Three column array | lon | lat | z |
% optional:
% fname - (String) Filename for .grd and .mask, no files will be
% created if no input
%
% OUTPUT
% lat - (1D Array) Latitude Array
% lon - (1D Array) Longitude Array
% z - (2D Array) Matrix of size (ny rows,nx columns)
%
% Created: John Goertz
% May 2013
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if xyz(2,1)-xyz(1,1) == 0
disp('Longitude and Latitude columns appear to be switched')
disp('Please make sure xyz columns are orderd | lon | lat | z | ')
return
else
dx = xyz(2,1) - xyz(1,1);
end
nx = round((max(xyz(:,1)) - min(xyz(:,1)))/dx + 1);
ny = round((max(xyz(:,2)) - min(xyz(:,2)))/dx + 1);
if max(xyz(:,1)) == xyz(nx,1)
else
disp('Check Longitude column in .xyz file for rounding errors')
return
end
lon = xyz(1:nx,1);
lat = xyz(1:nx:end,2);
z = zeros(ny,nx);
for ii = 1:ny
z(ii,(1:nx)) = xyz((1:nx)+(ii-1)*nx,3);
end
mask = ones(size(z));
loc = z >= 0;
mask(loc) = 0;
if exist('fname')
dlmwrite([fname,'.grd'],z,'precision','%.6f','delimiter', ' ')
dlmwrite([fname,'.mask'],mask,'delimiter',' ')
end
end