I'm having some trouble reading a WRF weather forecast prediction file where the coordinates are in a custom Lambert Conformal Conical projection.
My final goal is to get the weather forecast a a certain location given in a WGS84 projection (lat/lon), which I know for sure is within the NetCDF file.
I'm using a Python script and the NetCDF4 and PyProj libraries.
In the NetCDF file, downloaded from the local weather forecast service, there is the following information about the used projection:
<class 'netCDF4._netCDF4.Variable'> int32 Lambert_Conformal() grid_mapping_name: lambert_conformal_conic longitude_of_central_meridian: -14.100000381469727 standard_parallel: [43. 43.] latitude_of_projection_origin: 40.71900177001953 false_easting: -241.44951 false_northing: -4.10447 unlimited dimensions: current shape = () filling on, default _FillValue of -2147483647 used
So I try to recreate this using the pyproj library:
lat_0 = variables['Lambert_Conformal'].latitude_of_projection_origin lon_0 = variables['Lambert_Conformal'].longitude_of_central_meridian lat_1 = variables['Lambert_Conformal'].standard_parallel lat_2 = variables['Lambert_Conformal'].standard_parallel x_0 = variables['Lambert_Conformal'].false_easting y_0 = variables['Lambert_Conformal'].false_northing
I create a LCC origin projection with PyProj, as well as the destination WSG84 projection:
lcc = Proj(proj='lcc' , lat_1=lat_1, lat_2=lat_2, lat_0=lat_0, lon_0=lon_0, x_0=x_0, y_0=y_0, a=6370000, b=6370000) # Resulting projection string: '+proj=lcc +lat_0=40.7190017700195 +lon_0=-14.1000003814697 +lat_1=43 +lat_2=43 +x_0=-241.44951 +y_0=-4.10447 +R=6370000 +units=m +no_defs' wgs84 = Proj(proj='latlong', datum='WGS84')
Now that I have everything set-up, I perform the actual conversion between projections:
transformer = Transformer.from_proj(proj_from=wgs84, proj_to=lcc) # Geographical point which I know is contained within the weather forecast. lat_geo = 42.666981 lon_geo = -8.006540 coords=transformer.transform(lon_geo,lat_geo)
As result, I get:
which makes no sense to me, as I was expected something contained within the range of the variables in NetCDF file. For reference, the latitude variable in the NetCDF file is between 40.719 and 40.44141.
So, what am I missing here?
- Are the lat_0, lat_1, lat_2... parameters correctly loaded? Did I misplaced some of them?
- Is there some missing projection information in NetCDF file?
- I see that there is actually a single standard parallel (the same, but repeated twice). Does that make sense?
TLDR: How do I go from this custom LCC projection to standard WGS84 coordinates?
- The weather forecast agency says that the global lon/lat coordinates are supossed to be:
GeospatialCoverage: Longitude: -11.261331695726343 to -5.101928395629246 degrees_east Latitude: 40.421302033295966 to 45.230011083267385 degrees_north
- This is a map using their available web tools: