Re-projection Over the 180th Meridian
Anonym
With certain data formats, especially those with large swaths like AVHRR, NPP VIIIRS, and MODIS, there are times when the spatial extent of the data crosses over the 180th meridian. When measuring longitude in degrees for a projection, this can present a problem since the longitude is not continuous as it goes from -180 to 180 degrees. There are a few ways to correct for this anomaly in the data using ENVI and IDL in order to re-project it properly. For the examples below, let's say that we have the variables "lat","lon", and "data" read from an AVHRR HDF4 file. One way to get these variables using IDL:
e = envi(/current)
hdf_id = hdf_sd_start('my_hdf_file.hdf')
hdf_sd_fileinfo, hdf_id, numData,atts
dataset_names = strarr(numData)
for i=0, numData-1 do begin
hdf_sd_getinfo, hdf_sd_select(hdf_id, i),name=name
dataset_names[i] = name
endfor
index = where(dataset_names eq 'avhrr_band1')
lon_index = where(dataset_names eq 'longitude')
lat_index = where(dataset_names eq 'latitude')
dataset_id=hdf_sd_select(hdf_id, index)
hdf_sd_getdata, dataset_id, data
hdf_sd_endaccess, dataset_id
lon_id=hdf_sd_select(hdf_id, lon_index)
hdf_sd_getdata, lon_id, lon
hdf_sd_endaccess, lon_id
lat_id=hdf_sd_select(hdf_id, lat_index)
hdf_sd_getdata, lat_id, lat
hdf_sd_endaccess, lat_id
; Close the file:
hdf_sd_end,hdf_id
How can I tell if my data crosses the 180th meridan?
Depending on your dataset, there are few ways to test if you will need to correct for the extent crossing over from -180 degrees to +180 degrees. If you are working in a certain area, the US and Hawaii for example, you can test to see if your longitude is over a certain threshold. If you know that you are looking at data in North America, you can test to see if you have any longitude values near 180 degrees by taking the maximum of your longitude grid:
if max(lon)gt 179 then ...
A more robust way to check if your data set crosses the180th meridian is to take the difference between the minimum and maximum longitude. If the data does cross over the 180th meridian, the max will be ~180 and the minimum will be ~-180. It is unlikely that your data will wrap around the entire globe, like in the case of data around the North and South Poles, but that I'll save for another entry. So we can assume that if the difference between the max and min longitude is close to 360 degrees, the 180th meridian has been crossed:
if (max(lon)- min(lon)) gt 359 then ...
What can I do to correct for this case?
If you want to use the GLT reprojection tool available in ENVI to get the data on a regular grid, and your data crosses the 180th meridian, you are in luck! The engineers at here at Exelis put in the time to correct for this case automatically. If your data crosses the 180th meridian, the resulting image will be on the right side (+180 degrees) of a standard WGS-84 world map. If you'd rather it be on the left side of a WGS-84 map (-180 degrees), just subtract 360 degrees from the "lon" variable before doing GLT reprojection:
lon -= 360
Here is an example display in ENVI – in one image I applied the 360 shift. In the other I left it as is.
