Hi Tom,
Thanks for the background information.
> It should be possible to compute these vectors
It would be great if we could get compatibility across the netcdf-java and
proj4 stacks.
So, if the projection_x_coord and projection_y_coord units=m [ie "meters in a
plane tangent to the idealized earth nadir point"], the netcdf-java library
could then be asked to apply a transform to the coordinate data to view vectors
[radians] prior to reprojection calculations.
I can see that this can be done in the GEOSTransform.java; I'll have to read it
and the HRIT manual again to get a better idea of what's common and where the
transform to/from metres should go.
Cheers,
Leon
From: Tom Rink [mailto:rink@xxxxxxxxxxxxx]
Sent: Thursday, 9 July 2015 3:32 AM
To: netcdf-java@xxxxxxxxxxxxxxxx
Subject: Re: [netcdf-java] Geostationary CF grid_mapping problem: "HTTP Status
500 - Internal Server Error" [SEC=UNCLASSIFIED]
Hi Leon,
I worked with Martin, John Caron and others to define this standard, and
developed the software
on the Java side, so maybe I can help:
This transformation is formed by the cross product of two arithmetic
progressions, one in
W-E direction, the other S-N, forming a set of view vectors from an idealized
point in the
Earth's equatorial plane measured in radians. It should be possible to compute
these vectors
from a set of (x,y) measured in meters in a plane tangent to the idealized
earth nadir point
thus using most of the same software to support both transformation definitions.
Tom
On 7/8/15 3:03 AM, Leon Majewski wrote:
Hi Ryan, all, cc Martin,
I've found the issues, but they raise some questions.
How to get my gdal example to display WMS GetCapabilities:
According to cdm/src/main/java/ucar/nc2/constants/CF.java, the
grid_mapping_name attribute needs to be "geostationary", rather than
"geostationary_satellite" (as created by gdal)
ncatted -a
grid_mapping_name,geostationary_satellite,m,c,geostationary geos.nc
According to
cdm/src/main/java/ucar/unidata/geoloc/projection/sat/Geostationary.java
"fixed_angle_axis" or "sweep_angle_axis" is mandatory (not created by gdal)
ncatted -a sweep_angle_axis,geostationary_satellite,c,c,y
geos.nc
Need to rename so attributes to match code
ncrename -a
geostationary_satellite@satellite_height,perspective_point_height geos.nc
ncrename -a
geostationary_satellite@longitude_of_central_meridian,longitude_of_projection_origin
geos.nc
netcdf-java code doesn't seem to handle inverse flattening, change to
semi_minor_axis:
ncrename -a
geostationary_satellite@inverse_flattening,semi_minor_axis geos.nc
ncatted -a
semi_minor_axis,geostationary_satellite,m,f,6356752.3 geos.nc
Once these attributes are fixed, the layers are shown by the WMS
GetCapabilities request.
Great...
However, the GetMap doesn't work as expected
?FORMAT=image/jpeg&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&STYLES=boxfill/greyscale&WIDTH=800&HEIGHT=400&CRS=CRS:84&BBOX=0,-90,360,90&LAYERS=Band2&COLORSCALERANGE=0,255
Should effectively be the same as the result of:
gdalwarp -t_srs '+proj=latlong' -te 0 -90 360 90 -ts 800 400
NETCDF:geos.nc:Band1 out.tif -co "COMPRESS=JPEG" -ot Byte
(insert image... might not work)
[cid:image001.jpg@01D0BA1B.D0022830]
But the response only shows one pixel/colour.
This is because gdal, proj4 and pyproj... use units of 'm' for their
transforms, whereas this implementation seems to want to use radians (scan
angle) .
ncdump -v x geos.nc
x = -5498000 ... 5498000
Adding a guess at a scale factor for x, y
ncatted -a scale_factor,x,c,f,0.00000003 geos.nc
ncatted -a scale_factor,y,c,f,0.00000003 geos.nc
This time, the request shows some data, but obviously incorrect.
(insert image... might not work)
[cid:image002.jpg@01D0BA1B.D0022830]
So... there seems to be a difference in the units expected by GDAL/PROJ4 and
netcdf-java (and potentially CF). Is there a way to support "unit=m" within
netcdf-java and the CF standard? Alternatively proj4/gdal would need to be
updated to handle "units=rad"
I'm eager to hear your thoughts.
Cheers,
Leon
From: Leon Majewski
Sent: Tuesday, 19 May 2015 5:36 AM
To: 'Ryan May'
Cc: netcdf-java@xxxxxxxxxxxxxxxx<mailto:netcdf-java@xxxxxxxxxxxxxxxx>
Subject: RE: [netcdf-java] Geostationary CF grid_mapping problem: "HTTP Status
500 - Internal Server Error" [SEC=UNCLASSIFIED]
Hi Ryan,
I'll ask our admin to upgrade or spin up a vm. In the meantime, here's a simple
script to generate a file with similar issues (WMS page appears, but no layers).
# get some blue marble data (small):
wget
http://earthobservatory.nasa.gov/Features/BlueMarble/Images/land_shallow_topo_2048.tif
# add a projection
gdal_translate -a_srs '+proj=longlat' -a_ullr -180 90 180 -90
land_shallow_topo_2048.tif map.tif
# remap the data to our location (actually, use 80.5 due to dateline artefacts)
# output CF netcdf
gdalwarp -r cubic -tr 4000 4000 -t_srs "+proj=geos +lon_0=80.5 +h=35785863
+a=6378137.0 +b=6356752.3" -te -5500000 -5500000 5500000 5500000 map.tif
geos.nc -of NETCDF -co "COMPRESS=DEFLATE"
Cheers,
Leon
From: Ryan May [mailto:rmay@xxxxxxxx]<mailto:[mailto:rmay@xxxxxxxx]>
Sent: Tuesday, 19 May 2015 1:56 AM
To: Leon Majewski
Cc: netcdf-java@xxxxxxxxxxxxxxxx<mailto:netcdf-java@xxxxxxxxxxxxxxxx>
Subject: Re: [netcdf-java] Geostationary CF grid_mapping problem: "HTTP Status
500 - Internal Server Error" [SEC=UNCLASSIFIED]
Leon,
Version 4.6.1 of THREDDS was just released. Could you see if that fixes your
problem?
ftp://ftp.unidata.ucar.edu/pub/thredds/4.6/current/thredds.war
If not, it would be really helpful to get a sample of the problematic file so
we can try to reproduce locally.
Thanks,
Ryan
On Wed, May 13, 2015 at 12:06 AM, Leon Majewski
<Leon.Majewski@xxxxxxxxxx<mailto:Leon.Majewski@xxxxxxxxxx>> wrote:
Hi all,
> It's been suggested that I use the CF1.7 grid_mapping geostationary
I found that I need to set the grid_mapping_name to "geostationary_satellite"
float geostationary ;
geostationary:grid_mapping_name = "geostationary_satellite" ;
> HTTP Status 500
This fixes the HTTP Status 500. However, this results in their being no
queryable layers in the WMS and netcdf subset service fails [both worked with
vertical_perspective].
The WMS GetCapabilities for the vertical perspective version includes:
<Layer><Title>AGLS observations product suite</Title>
<Layer queryable="1"><Name>channel_0013_brightness_temperature</Name>
<Title>brightness temperature for channel 13 at 10.41 um</Title>
<Abstract>brightness temperature for channel 13 at 10.41 um</Abstract>
<EX_GeographicBoundingBox><westBoundLongitude>-180.0</westBoundLongitude><eastBoundLongitude>180.0</eastBoundLongitude><southBoundLatitude>-90.0</southBoundLatitude><northBoundLatitude>90.0</northBoundLatitude></EX_GeographicBoundingBox>
<BoundingBox CRS="CRS:84" minx="-180.0" maxx="180.0" miny="-90.0"
maxy="90.0"/>
<Dimension name="time" units="ISO8601" multipleValues="true"
current="true"
default="2015-02-20T00:00:20.000Z">2015-02-20T00:00:20.000Z</Dimension>
...
</Layer>
While the geostationary_satellite version has:
<Layer><Title>AGLS observations product suite</Title>
</Layer>
Has anyone see this behaviour/been able to work around it?
Thanks,
Leon
_______________________________________________
netcdf-java mailing list
netcdf-java@xxxxxxxxxxxxxxxx<mailto:netcdf-java@xxxxxxxxxxxxxxxx>
For list information or to unsubscribe, visit:
http://www.unidata.ucar.edu/mailing_lists/
--
Ryan May
Software Engineer
UCAR/Unidata
Boulder, CO
_______________________________________________
netcdf-java mailing list
netcdf-java@xxxxxxxxxxxxxxxx<mailto:netcdf-java@xxxxxxxxxxxxxxxx>
For list information or to unsubscribe, visit:
http://www.unidata.ucar.edu/mailing_lists/
