How to grid data with coordinates to create a spatial plot using geoshow

2017 8 31
2 个回答
回答已采纳
更新时间：2018 6 12
5 次查看（30 天）

0 个投票

I have 3 vectors: data, lat, lon which I am trying to plot spatially for the continental US. Is there a function where I can organize my lat and lon vectors in the appropriate gridded format which geoshow will plot properly while ensuring the data vector is organized in the same fashion so the data points remain with their respective coordinates?

请先登录，再回答此问题。

 采纳的回答

Chad Greene
Chad Greene 2017-8-31

1 个投票

Is it possible that your data vectors are in fact regular, but simply not in gridded format? To check, try
scatterm(lat,lon,20,data)
If the datapoints are regularly spaced, you should be able to use my xyz2grid function like this:
[LON,LAT,DATA] = xyz2grid(lon,lat,data);
where lowercase are your data vectors and upper case are gridded. Then you can do
pcolorm(LAT,LON,DATA)

13 个评论
显示 11更早的评论 隐藏 11更早的评论

unfortunately the data vectors are not regular; the scatterm function produced an error when running the program
Ronnie
Ronnie 2017-9-1
Error using gcm (line 25) Not a map axes.
Error in scatterm (line 48) gcm(ax);
Chad Greene
Chad Greene 2017-9-1
That means a map hasn't been initialized. Try
worldmap([min(lat) max(lat)],[min(lon) max(lon)])
scatterm(lat,lon,20,data)
Ronnie
Ronnie 2017-9-5
Thank you! This is working very well. I see the 20 denotes the area of each marker; do you know what units this is in? (i.e. if each point is symbolic of a 20 km x 20 km grid cell should I use 20?)
Chad Greene
Chad Greene 2017-9-5
The 20 isn't really units of real-world space. Rather, it's more like the number of pixels it takes up on your screen.
The scatterm suggestion was just to check whether your data correspond to a regular grid. Do the circles produced by scatterm make a gridded pattern? If yes, I recommend
% Grid the data:
[LON,LAT,DATA] = xyz2grid(lon,lat,data);
% Initialize a map:
worldmap([min(lat) max(lat)],[min(lon) max(lon)])
% plot the gridded data:
pcolorm(LAT,LON)
shading interp
Ronnie
Ronnie 2017-9-5
编辑：Walter Roberson 2017-9-6
I see what you are saying. Unfortunately my arrays are too large for the xyz2grid function to work. I am getting the following error message:
Error using accumarray
Requested 3857x5166498 (148.5GB) array exceeds maximum array size preference. Creation of
arrays greater than this limit may take a long time and cause MATLAB to become
unresponsive. See array size limit or preference panel for more information.
Error in xyz2grid (line 92)
Z = accumarray([yi xi],z(:),[],[],NaN);
Error in Monthly_Average_Function (line 85)
[X,Y,Z] = xyz2grid(double(Z),double(lon),double(lat));
I went into my preferences to uncheck the box that limits the array size limit to RAM; however once I did this MATLAB just shuts down. Do you have any suggestions on this issue?
Chad Greene
Chad Greene 2017-9-6
Oh wow, that's a lot of data points. I haven't downloaded the data you shared, but from KSSV's comment it looks like he found a lot of data points you can get rid of. If it's true, perhaps you can follow his suggestion
%%remove -9999
data(lon==-9999) = [] ;
lon(lon==-9999) = [] ;
lat(lat==-9999) = [] ;
And then use xyz2grid.
Ronnie
Ronnie 2017-9-6
编辑：Ronnie 2017-9-6
Even after removing these values I still get the same error message. Is there an advantage to using pcolor vs. scatterm ? I just want to ensure the spatial plot I am creating is an accurate representation
Walter Roberson
Walter Roberson 2017-9-6
编辑：Walter Roberson 2017-9-6
pcolor and scatterm are quite different purposes.
pcolor does a surface plot of a data grid -- it is surf() followed by view(2). pcolor of an N x M data grid produces (N-1) x (M-1) faces with the colors being interpolated by adjacent data points.
scatterm does a 2D scatter plot on a mapping axes, putting a marker at each location specified.
In answer to your earlier question, the area of the marker, such as 20, is in "points squared" where 1 point is 1/72 of an inch. In the case of circular markers, for size S, you would take R = sqrt(S/pi) to get the radius of the circle, in points. This size is used for the marker no matter how far in or out you zoom, so it does not reflect any particular data units.
Ronnie
Ronnie 2017-9-6
Thank you, Chad and Walter. I greatly appreciate your insights and I will move forward and utilize this knowledge.
Best, -Ronnie
Ronnie
Ronnie 2017-9-6
Update: I got it to work by trimming down the data even further (getting rid of NaN values) and then xyz2grid worked; and I am able to show the data using geoshow. Thank you so much for your help!

请先登录，再进行评论。

更多回答（1 个）

KSSV
KSSV 2017-8-31
编辑：Chad Greene 2017-8-31

1 个投票

Let data be your nx3 array which has lon, lat and data in the first, second and column respectively.
% Get longitude and latitude vectors
x = unique(data(:,1)) ;
y = unique(data(:,2)) ;
% dimensions of the data
nx = length(x) ;
ny = length(y) ;
% Frame matrix of grid
D = reshape(data(:,3),[ny,nx]) ;
% flip matrix to adjust for plot
H = flipud(H) ;
% Transpose the matrix
H = H' ; % Check if is required
surf(x,y,H) ;

16 个评论
显示 14更早的评论 隐藏 14更早的评论

In this solution, using the unique function to create latitude and longitude vectors will erase any duplicates of latitude and then any duplicates of longitude. It is okay to have multiple latitude readings, as long as they are paired with different respective longitude readings. I do not want to erase data points from the vectors. Both vecotrs (lat and lon) have already been updated to ensure that there are no duplicate coordinates.
I am seeking for a way to transform vectors in matrices which can be read by the function geoshow.
KSSV
KSSV 2017-8-31
The above code works..if your data is a structured grid and in (x,y,z) format......if it is unstructured grid, you need to follow the below:
x0 = min(lon) ; x1 = max(lon) ; nx = 100 ;
y0 = min(lat) ; y1 = max(lat) ; ny = 100 ;
x = linspace(x0,x1,nx) ;
y = linspace(y0,y1,ny) ;
[X,Y] = meshgrid(x,y) ;
Z = griddata(lon,lat,data,X,Y)
surf(X,Y,Z)
Ronnie
Ronnie 2017-8-31
My data is not a structured grid. And I cannot use meshgrid for my lat and lon vectors unfortunately because the vector sizes I am using gives me an error saying it is too large to use in meshgrid. I appreciate your help!
Ronnie
Ronnie 2017-8-31
I am also looking to do a 2d plot, thus the surf function will not work, which is why I am trying to use geoshow.
I appreciate your help. I will send the data tomorrow once I am back at the computer!
Ronnie
Ronnie 2017-9-1
SMAP_March was the 3-columned matrix with data in column 1 and the data's corresponding lat and lon in columns 2 and 3 respectively.
data=SMAP_March(:,1); lat=SMAP_March(:,2); lon=SMAP_March(:,3);
%Saves a variable as the script file (will replace rest of script file) matlab.io.saveVariablesToScript('DATA.m','data') matlab.io.saveVariablesToScript('LAT.m','lat') matlab.io.saveVariablesToScript('LON.m','lon')
I used the following code to save the 3 vectors I am sending to you (attached).
KSSV
KSSV 2017-9-2
Week end...I will get back to you on Monday ..:)
Enjoy the weekend, I will not have access to my computer until Monday as well. :)
The site will not allow me to send any files greater than 5 MB, which the lat and data files are. The LON.mat file was bellow the size limit and is attached. Can you please send me an email at ronnie.aggie2016@gmail.com and I will respond with the .mat files. I apologize for the inconvenience
KSSV
KSSV 2017-9-5
You may attach the file in your google drive and past the link here......
KSSV
KSSV 2017-9-6
This works....
lon = load('LON.mat') ;
lon = lon.lon ;
lat = load('LAT.mat') ;
lat = lat.lat ;
data = load('DATA.mat') ;
data = data.data ;
%%remove -9999
data(lon==-9999) = [] ;
lon(lon==-9999) = [] ;
lat(lat==-9999) = [] ;
scatter(lon,lat,data,data)
Ronnie
Ronnie 2017-9-6
Thank you for the response! Much appreciation.
-Ronnie
TAPAS
TAPAS 2018-6-12
The code xyz2grid is not working it's showing mistake in line 31 in xyz read and line 72 in in xyz2grid

请先登录，再进行评论。

类别

在 帮助中心 和 File Exchange 中查找有关 Lengths and Angles 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by