File Exchange

image thumbnail

Find pixel indices in HDF-EOS files based on LatLon coordinates

version 1.1.0.0 (5.8 KB) by M Sohrabinia
Finds x-y direct indexes of one or more LatLon coordinate pairs in certain HDF-EOS files

5 Downloads

Updated 15 Jun 2012

View License

This function takes Latitude-Longitude pairs with name of an HDF-EOS file and finds the X-Y indexes of each LatLon pair in the HDF file. The function will work with HDF files which have:
1. Sinusoidal Projection (e.g., MODIS LST Level 3 product);
2. Global Cylindrical Equal-Area Projection or Cylindrical Equal-Area Scalable Earth Grid (EASE)
Important: LatLon pairs must be given in a matrix where first column contains Latitude and the SECOND COLUMN LONGITUDE, while the function will write X as first and Y as the second output. One or more pairs of LatLon coordinates can be given, where each line will contain coordinates of one point. The first two output variables will contain the X and Y indexes with the same order as LatLon inputs. The function is tested with the MODIS LST Level 3 (MOD11A1) and MODIS BRDF/Albedo (MCD43B3) products for proj=1 (Sinusoidal), and AMSR-E Soil Moisture (AMSR_E_L3_DailyLand_V07) product for proj=2 (EASE grid).
--Input variables:
LatLon: matrix with two columns: Latitudes and Longitudes
filename: HDF file name (including filepath if file is not in current
directory)
proj: Projection type of the HDF file (currently only 1 or 2)

--Output variables:
coordX_idx: x indexes of given LatLon coordinate(s) in the given HDF file
coordY_idx: y indexes of given LatLon coordinate(s) in the given HDF file (if
0 is returned by coordX_idx or coordY_idx means the LatLon pair could
not be found in the given image/HDF file).
projType: (string) full name of the given projection type
mdist2point: minimum distance to each point in meters that the function
could work out (should not be bigger than xRes or yRes).
xRes: resolution of the image along x or longitudes (it is always positive)
yRes: resolution of the image along y or latitudes (since Latitudes decrease from UppeLeft to LowerRight corner of an image, this will be a negative value--> note: this is now updated to be always positive).

--Example:
LatLon=[-43.5,172.5; -43.9,172.75]
%download the following file and put it in the current directory:
ftp://e4ftl01.cr.usgs.gov/MODIS_Dailies_E/MOLA/MYD11A1.005/2012.05.07/MYD11A1.A2012128.h30v13.005.2012129204401.hdf
filename='./MYD11A1.A2012128.h30v13.005.2012129204401.hdf';
proj=1;
[coordX_idx,coordY_idx,projType,mdist2point,xRes,yRes] = ...
findIndexes(LatLon,filename, proj)

The returned coordX_idx coordY_idx are indexes of the LST pixels covering 1km^2 in the above coordinates.

First Version: Feb 19 2012 (V01)
Added functionailty: Apr 25 2012 (V02)
Updated: Jun 15 2012

Comments and Ratings (3)

ICY LEE

Thank you for your code.
I have a question that the coordX_idx,coordY_idx that reprsent the row and colomn of the data in matlab?I found X represent colomn and Y represent row according to ARCGIS.
What method do you use to ensure the result?In Arcgis the point of the same latlon don't always have the same value.

Shanshan

Line 202 and 203 do not work for MODIS tiles with sinusoidal projection in North America. An option is to change the two lines as:
x=abs(x*dist1degree-ULx0);
y=abs(y*dist1degree-ULy0);

Did not test for other places or projections.

Shanshan

Updates

1.1.0.0

1. Coordinates of the pixels to find are now matched with the center of each pixel in contrast to the UpperLeft corner in the previous code.
2. yRes is changed to be always positive, description of the function is updated.

MATLAB Release Compatibility
Created with R2008b
Compatible with any release
Platform Compatibility
Windows macOS Linux