[ncl-talk] Regridding from OSGB spatial coordinate to Latitude-longitude in rotated pole coordinates

Dennis Shea shea at ucar.edu
Thu Feb 25 14:07:35 MST 2021


I don't understand:

"This data *can't also read* correctly with Climate Data Operator (CDO)"

To me something is not correct with the file if the CDO and NCL are having
issues.
Python is rendering a correct plot?    Somebody has programmed around the
issue.
===

You can go back and forth between the lat/lon gridded coordinates and x/y
as indicated by the equations at:

*https://idlecoding.com/transverse-mercator-with-python/comment-page-1/#:~:text=%20Transverse%20Mercator%20With%20Python%20%201%20The,any%20code%20to%20transform%20an%20actual...%20More%20*
<https://idlecoding.com/transverse-mercator-with-python/comment-page-1/#:~:text=%20Transverse%20Mercator%20With%20Python%20%201%20The,any%20code%20to%20transform%20an%20actual...%20More%20>

But to me, these are both on the file.

        double grid_longitude(projection_y_coordinate,
projection_x_coordinate) ;
        double grid_latitude(projection_y_coordinate,
projection_x_coordinate) ;

        double projection_x_coordinate(projection_x_coordinate) ;
        double projection_y_coordinate(projection_y_coordinate) ;

The information at:

       int transverse_mercator ;
                transverse_mercator:grid_mapping_name =
"transverse_mercator" ;
                transverse_mercator:longitude_of_prime_meridian = 0. ;
                transverse_mercator:semi_major_axis = 6377563.396 ;
                transverse_mercator:semi_minor_axis = 6356256.909 ;
                transverse_mercator:longitude_of_central_meridian = -2. ;
                transverse_mercator:latitude_of_projection_origin = 49. ;
                transverse_mercator:false_easting = 400000. ;
                transverse_mercator:false_northing = -100000. ;
                transverse_mercator:scale_factor_at_central_meridian =
0.9996012717 ;

is possibly useful, I have tried to use two of the attributes but it did
not work.




On Thu, Feb 25, 2021 at 1:08 PM S Br <sbr.climate at gmail.com> wrote:

> Hi Dennis,
> The plot is perfect with Python. Please see the attached image.
> The data I am talking about is generated using Python.
>
> This data can't also read correctly with Climate Data Operator (CDO)
>
> Thanks.
>
>
> On Thu, Feb 25, 2021 at 7:26 PM Dennis Shea <shea at ucar.edu> wrote:
>
>> Obviously, not over the UK but ......
>>
>> I think more information from you might be necessary.
>> "transverse_mercator" is not supported by NCL options.
>>
>> You state: "The data is correctly displayed in Python but not in NCL."
>>
>> Well, what does the plot look like via Python?
>>
>> On Thu, Feb 25, 2021 at 12:19 PM Dennis Shea <shea at ucar.edu> wrote:
>>
>>> I **speculate** that the attributes associated with the
>>> "transverse_mercator" variable must be used.
>>> ====================================================
>>> ;---Get the source lat/lon grid
>>> tas   = f->tasmax(0,0,:,:)
>>>
>>> lat2d = f->grid_latitude
>>> lon2d = f->grid_longitude
>>> printMinMax(lat2d,0)
>>> printMinMax(lon2d,0)
>>>
>>> itm   = f->transverse_mercator
>>>
>>> print(itm)
>>> print("itm at latitude_of_projection_origin
>>> ="+itm at latitude_of_projection_origin)
>>> print("itm at longitude_of_central_meridian
>>> ="+itm at longitude_of_central_meridian)
>>>
>>> ; clarity
>>>
>>> lat_shift = itm at latitude_of_projection_origin
>>> lon_shift = itm at longitude_of_central_meridian
>>>
>>> ; apply 'shift'
>>>
>>> lat2d = lat2d + lat_shift
>>> lon2d = lon2d + lon_shift
>>> printMinMax(lat2d,0)
>>> printMinMax(lon2d,0)
>>>
>>> tas at lat2d = lat2d             ; use 'shifted' coordinates
>>> tas at lon2d = lon2d
>>> ================
>>>
>>> See attached
>>>
>>>
>>> On Thu, Feb 25, 2021 at 12:00 PM Dave Allured - NOAA Affiliate via
>>> ncl-talk <ncl-talk at mailman.ucar.edu> wrote:
>>>
>>>> I see what's going on.  From your other e-mail, the data file contains
>>>> two different sets of coordinates.  There are the projection coordinates in
>>>> meters.  Those are not latitudes and longitudes, therefore not sufficient
>>>> for NCL.  But the geographic coordinates are also provided.  Double
>>>> coordinates are frequently encountered in GIS related applications.  Try
>>>> this:
>>>>
>>>>     tas at lat2d=f->grid_latitude
>>>>     tas at lon2d=f->grid_longitude
>>>>
>>>>
>>>> On Thu, Feb 25, 2021 at 11:42 AM S Br <sbr.climate at gmail.com> wrote:
>>>>
>>>>> If I print tas at lat2d and tas at lon2d, they give the following
>>>>> information.
>>>>>
>>>>> Variable: taslat
>>>>> Type: double
>>>>> Total Size: 896 bytes
>>>>>             112 values
>>>>> Number of Dimensions: 1
>>>>> Dimensions and sizes: [112]
>>>>> Coordinates:
>>>>> (0) -102000
>>>>> (1) -90000
>>>>> (2) -78000
>>>>> (3) -66000
>>>>> (4) -54000
>>>>> (5) -42000
>>>>> (6) -30000
>>>>> (7) -18000
>>>>> (8) -6000
>>>>> (9) 6000
>>>>> (10) 18000
>>>>> (11) 30000
>>>>> (12) 42000
>>>>> (13) 54000
>>>>> (14) 66000
>>>>> (15) 78000
>>>>> (16) 90000
>>>>> (17) 102000
>>>>> (18) 114000
>>>>> (19) 126000
>>>>> (20) 138000
>>>>> (21) 150000
>>>>> (22) 162000
>>>>> (23) 174000
>>>>> (24) 186000
>>>>> (25) 198000
>>>>> (26) 210000
>>>>> (27) 222000
>>>>> (28) 234000
>>>>> (29) 246000
>>>>> (30) 258000
>>>>> (31) 270000
>>>>> (32) 282000
>>>>> (33) 294000
>>>>> (34) 306000
>>>>> (35) 318000
>>>>> (36) 330000
>>>>> (37) 342000
>>>>> (38) 354000
>>>>> (39) 366000
>>>>> (40) 378000
>>>>> (41) 390000
>>>>> (42) 402000
>>>>> (43) 414000
>>>>> (44) 426000
>>>>> (45) 438000
>>>>> (46) 450000
>>>>> (47) 462000
>>>>> (48) 474000
>>>>> (49) 486000
>>>>> (50) 498000
>>>>> (51) 510000
>>>>> (52) 522000
>>>>> (53) 534000
>>>>> (54) 546000
>>>>> (55) 558000
>>>>> (56) 570000
>>>>> (57) 582000
>>>>> (58) 594000
>>>>> (59) 606000
>>>>> (60) 618000
>>>>> (61) 630000
>>>>> (62) 642000
>>>>> (63) 654000
>>>>> (64) 666000
>>>>> (65) 678000
>>>>> (66) 690000
>>>>> (67) 702000
>>>>> (68) 714000
>>>>> (69) 726000
>>>>> (70) 738000
>>>>> (71) 750000
>>>>> (72) 762000
>>>>> (73) 774000
>>>>> (74) 786000
>>>>> (75) 798000
>>>>> (76) 810000
>>>>> (77) 822000
>>>>> (78) 834000
>>>>> (79) 846000
>>>>> (80) 858000
>>>>> (81) 870000
>>>>> (82) 882000
>>>>> (83) 894000
>>>>> (84) 906000
>>>>> (85) 918000
>>>>> (86) 930000
>>>>> (87) 942000
>>>>> (88) 954000
>>>>> (89) 966000
>>>>> (90) 978000
>>>>> (91) 990000
>>>>> (92) 1002000
>>>>> (93) 1014000
>>>>> (94) 1026000
>>>>> (95) 1038000
>>>>> (96) 1050000
>>>>> (97) 1062000
>>>>> (98) 1074000
>>>>> (99) 1086000
>>>>> (100) 1098000
>>>>> (101) 1110000
>>>>> (102) 1122000
>>>>> (103) 1134000
>>>>> (104) 1146000
>>>>> (105) 1158000
>>>>> (106) 1170000
>>>>> (107) 1182000
>>>>> (108) 1194000
>>>>> (109) 1206000
>>>>> (110) 1218000
>>>>> (111) 1230000
>>>>>
>>>>>
>>>>> Variable: taslon
>>>>> Type: double
>>>>> Total Size: 656 bytes
>>>>>             82 values
>>>>> Number of Dimensions: 1
>>>>> Dimensions and sizes: [82]
>>>>> Coordinates:
>>>>> (0) -210000
>>>>> (1) -198000
>>>>> (2) -186000
>>>>> (3) -174000
>>>>> (4) -162000
>>>>> (5) -150000
>>>>> (6) -138000
>>>>> (7) -126000
>>>>> (8) -114000
>>>>> (9) -102000
>>>>> (10) -90000
>>>>> (11) -78000
>>>>> (12) -66000
>>>>> (13) -54000
>>>>> (14) -42000
>>>>> (15) -30000
>>>>> (16) -18000
>>>>> (17) -6000
>>>>> (18) 6000
>>>>> (19) 18000
>>>>> (20) 30000
>>>>> (21) 42000
>>>>> (22) 54000
>>>>> (23) 66000
>>>>> (24) 78000
>>>>> (25) 90000
>>>>> (26) 102000
>>>>> (27) 114000
>>>>> (28) 126000
>>>>> (29) 138000
>>>>> (30) 150000
>>>>> (31) 162000
>>>>> (32) 174000
>>>>> (33) 186000
>>>>> (34) 198000
>>>>> (35) 210000
>>>>> (36) 222000
>>>>> (37) 234000
>>>>> (38) 246000
>>>>> (39) 258000
>>>>> (40) 270000
>>>>> (41) 282000
>>>>> (42) 294000
>>>>> (43) 306000
>>>>> (44) 318000
>>>>> (45) 330000
>>>>> (46) 342000
>>>>> (47) 354000
>>>>> (48) 366000
>>>>> (49) 378000
>>>>> (50) 390000
>>>>> (51) 402000
>>>>> (52) 414000
>>>>> (53) 426000
>>>>> (54) 438000
>>>>> (55) 450000
>>>>> (56) 462000
>>>>> (57) 474000
>>>>> (58) 486000
>>>>> (59) 498000
>>>>> (60) 510000
>>>>> (61) 522000
>>>>> (62) 534000
>>>>> (63) 546000
>>>>> (64) 558000
>>>>> (65) 570000
>>>>> (66) 582000
>>>>> (67) 594000
>>>>> (68) 606000
>>>>> (69) 618000
>>>>> (70) 630000
>>>>> (71) 642000
>>>>> (72) 654000
>>>>> (73) 666000
>>>>> (74) 678000
>>>>> (75) 690000
>>>>> (76) 702000
>>>>> (77) 714000
>>>>> (78) 726000
>>>>> (79) 738000
>>>>> (80) 750000
>>>>> (81) 762000
>>>>>
>>>>> On Thu, Feb 25, 2021 at 6:26 PM Dave Allured - NOAA Affiliate <
>>>>> dave.allured at noaa.gov> wrote:
>>>>>
>>>>>> It is starting to look like there is an error in the attached
>>>>>> coordinates.  What are the min and max values of tas at lat2d and
>>>>>> tas at lon2d?
>>>>>>
>>>>>>
>>>>>> On Thu, Feb 25, 2021 at 10:45 AM S Br <sbr.climate at gmail.com> wrote:
>>>>>>
>>>>>>> Hi Dennis,
>>>>>>> The data draws on a wrong location and the map boundary is not
>>>>>>> aligning.
>>>>>>> It should plot over the UK which is about 50-60N.
>>>>>>>
>>>>>>> Thanks.
>>>>>>>
>>>>>>> On Thu, Feb 25, 2021 at 5:41 PM Dennis Shea <shea at ucar.edu> wrote:
>>>>>>>
>>>>>>>> It helps to give NCL the correct lat/lon array information.
>>>>>>>>
>>>>>>>> See attached
>>>>>>>>
>>>>>>>> On Thu, Feb 25, 2021 at 10:23 AM Dave Allured - NOAA Affiliate via
>>>>>>>> ncl-talk <ncl-talk at mailman.ucar.edu> wrote:
>>>>>>>>
>>>>>>>>> That script looks pretty good.  You can see from the error
>>>>>>>>> messages that NCL is trying to understand @lat2d and @lon2d, but something
>>>>>>>>> is not quite right with them.  We need to see printVarSummary for tas,
>>>>>>>>> tas at lat2d, and tas at lon2d.  It is possible that dimensions are
>>>>>>>>> mismatched between data and coordinates, for some reason.
>>>>>>>>>
>>>>>>>>> There is also this note in that doc page:  "If your data variable
>>>>>>>>> is a two-dimensional (2D) array ordered lon x lat, you will need to reorder
>>>>>>>>> it to be lat x lon."  Please check and see if this is the problem.
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> On Thu, Feb 25, 2021 at 10:00 AM S Br <sbr.climate at gmail.com>
>>>>>>>>> wrote:
>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> Hi Dave,
>>>>>>>>>> Many thanks for your reply. Sorry I missed to provide the details.
>>>>>>>>>> I have run the following NCL script  and came across the errors
>>>>>>>>>> given below.
>>>>>>>>>> You can see the attached plot that it draws value throughout the
>>>>>>>>>> globe. It is supposed to be drawn only over the UK.
>>>>>>>>>> Am I doing something wrong?
>>>>>>>>>>
>>>>>>>>>> ;********
>>>>>>>>>> srcFileName="tas_rcp85_land-rcm_uk_12km_01_mon_198012-208011.nc"
>>>>>>>>>>
>>>>>>>>>> ;---Read File
>>>>>>>>>> f=addfile(srcFileName,"r")
>>>>>>>>>> print(f)
>>>>>>>>>>
>>>>>>>>>> ;---Get the source lat/lon grid
>>>>>>>>>> tas=f->tas(0,0,:,:)
>>>>>>>>>> tas at lat2d=f->projection_y_coordinate
>>>>>>>>>> tas at lon2d=f->projection_x_coordinate
>>>>>>>>>>
>>>>>>>>>> ;---Plot data
>>>>>>>>>>   wks = gsn_open_wks("x11","map")
>>>>>>>>>>
>>>>>>>>>>   res              = True
>>>>>>>>>>   res at gsnMaximize  = True     ; maximize plot in frame
>>>>>>>>>>   res at cnFillOn     = True     ; turn on contour fill
>>>>>>>>>>   res at cnLinesOn    = False    ; turn off contour fill
>>>>>>>>>>   res at gsnAddCyclic     = False
>>>>>>>>>> ;   res at mpMinLatF    =  48.0             ; only plot 30S to 30N
>>>>>>>>>> ;   res at mpMaxLatF    =  60.0
>>>>>>>>>> ;  res at mpMinLonF    =  -15.0          ; only plot 30S to 30N
>>>>>>>>>> ;   res at mpMaxLonF    = 6.0
>>>>>>>>>>  res at mpDataBaseVersion = "MediumRes"
>>>>>>>>>>    res at mpLandFillColor = "white"
>>>>>>>>>>
>>>>>>>>>>   plot = gsn_csm_contour_map(wks,tas(:,:),res)
>>>>>>>>>>
>>>>>>>>>> ********************
>>>>>>>>>> (0) is_valid_latlon2d_attr: Warning: The 'lat2d' attribute must
>>>>>>>>>> either be
>>>>>>>>>> (0) the same dimension sizes as the data, or one element larger
>>>>>>>>>> in both directions.
>>>>>>>>>> (0) Your data will most likely not be overlaid on the map
>>>>>>>>>> correctly.
>>>>>>>>>> (0) check_for_y_lat_coord: Warning: Data either does not contain
>>>>>>>>>> (0) a valid latitude coordinate array or doesn't contain one at
>>>>>>>>>> all.
>>>>>>>>>> (0) A valid latitude coordinate array should have a 'units'
>>>>>>>>>> (0) attribute equal to one of the following values:
>>>>>>>>>> (0)    'degrees_north' 'degrees-north' 'degree_north' 'degrees
>>>>>>>>>> north' 'degrees_N' 'Degrees_north' 'degree_N' 'degreeN' 'degreesN' 'deg
>>>>>>>>>> north'
>>>>>>>>>> (0) is_valid_latlon2d_attr: Warning: The 'lon2d' attribute must
>>>>>>>>>> either be
>>>>>>>>>> (0) the same dimension sizes as the data, or one element larger
>>>>>>>>>> in both directions.
>>>>>>>>>> (0) Your data will most likely not be overlaid on the map
>>>>>>>>>> correctly.
>>>>>>>>>> (0) check_for_lon_coord: Warning: Data either does not contain
>>>>>>>>>> (0) a valid longitude coordinate array or doesn't contain one at
>>>>>>>>>> all.
>>>>>>>>>> (0) A valid longitude coordinate array should have a 'units'
>>>>>>>>>> (0) attribute equal to one of the following values:
>>>>>>>>>> (0)    'degrees_east' 'degrees-east' 'degree_east' 'degrees east'
>>>>>>>>>> 'degrees_E' 'Degrees_east' 'degree_E' 'degreeE' 'degreesE' 'deg east'
>>>>>>>>>>
>>>>>>>>>> On Thu, Feb 25, 2021 at 4:17 PM Dave Allured - NOAA Affiliate <
>>>>>>>>>> dave.allured at noaa.gov> wrote:
>>>>>>>>>>
>>>>>>>>>>> Does your file contain its own lat and lon coordinate arrays?
>>>>>>>>>>> You will need to check this yourself with ncdump or NCL plus print commands.
>>>>>>>>>>>
>>>>>>>>>>> If the file has its own coordinate arrays, then do not regrid or
>>>>>>>>>>> reproject.  Doing that will lose detail and reduce the quality of the
>>>>>>>>>>> plotted image.  NCL is quite capable of directly displaying any data with
>>>>>>>>>>> attached coordinate arrays, and it does not need to know anything about the
>>>>>>>>>>> projection of the data.  Please see documentation for "Plotting data on a
>>>>>>>>>>> map" for much more information about aligning data on a map plot.
>>>>>>>>>>>
>>>>>>>>>>>     https://www.ncl.ucar.edu/Applications/plot_data_on_map.shtml
>>>>>>>>>>>
>>>>>>>>>>> If the coordinate arrays in the file are two dimensional, then
>>>>>>>>>>> see example 3, "Curvilinear grid".  When you have 2-D coordinates, they
>>>>>>>>>>> must be passed to the NCL plot routine with special attributes
>>>>>>>>>>> *lat2d* and *lon2d*, or alternatively *sfXArray* and *sfYArray* attached
>>>>>>>>>>> to the plot resource variable.  Note that 1-D coordinates do not need this
>>>>>>>>>>> special treatment because they are normally attached directly to the data
>>>>>>>>>>> variable to plot.
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> On Thu, Feb 25, 2021 at 8:24 AM S Br via ncl-talk <
>>>>>>>>>>> ncl-talk at mailman.ucar.edu> wrote:
>>>>>>>>>>>
>>>>>>>>>>>> Hi All,
>>>>>>>>>>>> I have a set of regional gridded data in the British National
>>>>>>>>>>>> Grid (OSGB) spatial coordinate system. I am looking to use this data in NCL
>>>>>>>>>>>> but the projection is not correctly displayed? Could it be possible to
>>>>>>>>>>>> regrid this data to Latitude-longitude in rotated pole coordinates or to a
>>>>>>>>>>>> regular lat/lon grid.
>>>>>>>>>>>>
>>>>>>>>>>>> I would like to mention that these data sets are generated
>>>>>>>>>>>> using the Iris package in Python. The data is correctly displayed in Python
>>>>>>>>>>>> but not in NCL.
>>>>>>>>>>>>
>>>>>>>>>>>> For your convenience I have attached here a sample NetCDF file.
>>>>>>>>>>>>
>>>>>>>>>>>> Thank you.
>>>>>>>>>>>> SBR
>>>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>> ncl-talk mailing list
>>>> ncl-talk at mailman.ucar.edu
>>>> List instructions, subscriber options, unsubscribe:
>>>> https://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>
>>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20210225/21979327/attachment.html>


More information about the ncl-talk mailing list