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

S Br sbr.climate at gmail.com
Fri Feb 26 02:43:44 MST 2021


Hi Karin,
Thanks for your comment. Yes, I had posted this problem also in the CDO
forum.
Hope you can find some solution using CDO.

Thanks.

On Fri, Feb 26, 2021 at 8:42 AM Karin Meier-Fleischer via ncl-talk <
ncl-talk at mailman.ucar.edu> wrote:

> Hi,
>
> here, I weigh in to make some comments about CDO and the Python Iris
> projections. Your question about the OSGB data in the CDO forum had me
> running into the same problems days ago.
>
> CDO usually can handle pre-projected grids using a gridfile with the
> grid_mapping informaion, but in your case there is something missing to
> create a correct gridfile. Transverse Mercator projection is quite special.
>
> The python package Iris provides an extra class for the projection type
> OSGB.
>
>    class iris.coord_systems.OSGBA
>       Specific transverse mercator projection on a specific ellipsoid
>
> I'll see how Iris has solved the problem and get back to you, but maybe
> not today;).
>
> Ciao,
> Karin
>
> Am 25.02.21 um 22:07 schrieb Dennis Shea via ncl-talk:
> > 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
> > <mailto: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
> >     <mailto: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
> >         <mailto: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
> >             <mailto: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 <mailto: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
> >                     <mailto: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
> >                         <mailto: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 <mailto: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
> >                                 <mailto: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
> >                                     <mailto: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
> >                                         <
> http://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
> >                                         <mailto: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
> >                                             <
> 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
> >                                             <mailto:
> 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 <mailto:
> ncl-talk at mailman.ucar.edu>
> >                 List instructions, subscriber options, unsubscribe:
> >                 https://mailman.ucar.edu/mailman/listinfo/ncl-talk
> >                 <https://mailman.ucar.edu/mailman/listinfo/ncl-talk>
> >
> >
> > _______________________________________________
> > ncl-talk mailing list
> > ncl-talk at mailman.ucar.edu
> > List instructions, subscriber options, unsubscribe:
> > https://mailman.ucar.edu/mailman/listinfo/ncl-talk
> >
> _______________________________________________
> 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/20210226/9a346077/attachment.html>


More information about the ncl-talk mailing list