How to read ASCII grid format

20 views (last 30 days)
Melissa
Melissa on 13 Jan 2015
Commented: Stephen23 on 24 Feb 2016
Hello,
I am confused about how to read this ASCII grid I have.
Let's call it A.asc A has a resolution of 0.5 x 0.5 degrees
ncols 720
nrows 360
xllcorner -180.25
yllcorner -89.75
cellsize 0.5
NODATA_value -9999
I am confused about the gridding format because this dataset should be worldwide and have even intervals of 0.5
How should I make the latitude and longitude vectors?
If I use linspace, I will not get a 360x1 vector evenly spaced by 0.5
Does it go from -180.25 to 0? If so, then what happened to the coordinates between 0 and 180?
Thank you,
Melissa
  3 Comments
Melissa
Melissa on 13 Jan 2015
Ah, yes. I, for some reason, was relating the longitude with the 360 dimension haha.
Yes, the file is (360x720 double)
However, if I try to make a longitude vector with linspace
A = linspace(-180.25,179.75,360);
Then I get a 1x360 dimension vector, but the numbers within them will not be spaced exactly 0.5 degrees apart from each other. Is there another way to do this?
Stephen23
Stephen23 on 24 Feb 2016
Bouchra AITHSSAINE's "Answer" moved here:
Hi , i wanna ask you the same question , I want to extract a matrix based on a large matrix that I have in ascii format, the problem is that I dont have the number of row or column of the matrix I just have the geographic coordinates, so how to do? Thank you

Sign in to comment.

Accepted Answer

Mohammad Abouali
Mohammad Abouali on 13 Jan 2015
Edited: Mohammad Abouali on 14 Jan 2015
Since you already have Mapping Toolbox, the easiest is to use pix2latlon(). Something like this:
% Read your data
[Data,RefMat]= arcgridread('youtfile.shp');
% prepare a lat/lon grid
[nrows,ncols,~]=size(Data);
[row,col]=ndgrid(1:nrows,1:ncols);
[lat,lon]=pix2latlon(RefMat,row,col);
% NOTE: lat(:,1) and lon(1,:)' are the vectors that you are looking
% display a sample field
figure,
worldmap([min(lat(:)) max(lat(:))],[min(lon(:)) max(lon(:))])
geoshow(lat,lon,Data(:,:,1),'DisplayType','surface');
% display on a globe (more fancy)
figure
axesm('globe','Geoid',earthRadius)
meshm(Data(:,:,1), RefMat)
view(3)
set(gca,'Box','off')
  1 Comment
Mohammad Abouali
Mohammad Abouali on 14 Jan 2015
Edited: Mohammad Abouali on 14 Jan 2015
% NOTE: lat(:,1) and lon(1,:)' are the vectors that you are looking

Sign in to comment.

More Answers (1)

Mohammad Abouali
Mohammad Abouali on 13 Jan 2015
Edited: Mohammad Abouali on 13 Jan 2015
A bit of explanation for ESRI Arc ASCII Grid Format:
At the end you get an image of with nRows x nCols So it means your data is 360x720.
If you check after No Data there should be n*360 lines and 720 number on each line. (n is number of bands).
xll,yll, and cell size are obvious.
I usually turn No-Data to NaN. NODATA_value means anywhere you see that number it wasn't any data for it.
Note that the ascii file is written like an image. So the left most number on the line is for on the west side and the right most number on the line is for the east side.
Also first line is on top (so north) and the last line in the file is located on the south part of your image. Therefore your lat/lon or coordinates would be
lon=(0:-1:-(nrows-1)).*cellsize + yllcorner + cellsize/2
lat=(0:ncols).*celsize + xllcorner + cellsize/2
note the cellsize/2 because x/yllcorner is the lower left corner of the cell. If instead of x/yllcorner you had x/yllcenter then it was the center of the cell.
You longitude seems to be from -180 to 179.5 (cell centers) So the entire globe.
But I see that your yllcorner is -89.75. That is kinda strange. If it is indeed like that some how your ascii file is written from south to north, which is not usual. Check if it was -89.75 or +89.75.
You can write your own reader, it is very easy. Or you can use the MATLAB's arcgridread().
  3 Comments
Mohammad Abouali
Mohammad Abouali on 13 Jan 2015
Edited: Mohammad Abouali on 13 Jan 2015
That's the reference matrix. if you have the row and col number of a pixel then the lat/lon coordinate would be:
[lon lat] = [row col 1] * R
You can use it to calculate the coordinate of each pixel, or show the data on the coordinate. Check mapshow() or geoshow() for some examples.
That's strange though that the yllcorner starts from the south. ESRI ASCII standard says that it should start from north. Although, I have to say that Lat=90.5 doesn't exist.
If your original question is answered please accept the answer. If you have further question, please start a new questions.
Melissa
Melissa on 13 Jan 2015
Edited: Melissa on 13 Jan 2015
The geomap command is very helpful. However, I would like to be able to create a vector for both the latitude and longitude.
Latitude = 180x1 double
Longitude = 360x1 double
So that I could perform operations with the data attached to them. How can I accomplish this with the reference matrix given by arcgridread()? I'm not sure what you mean about the row and column number of a pixel

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!