1212from mpas_tools .mesh .creation .jigsaw_driver import jigsaw_driver
1313from mpas_tools .mesh .creation .jigsaw_to_netcdf import jigsaw_to_netcdf
1414from mpas_tools .viz .colormaps import register_sci_viz_colormaps
15+ from mpas_tools .logging import LoggingContext
1516
1617
1718def build_spherical_mesh (cellWidth , lon , lat , earth_radius ,
1819 out_filename = 'base_mesh.nc' , plot_cellWidth = True ,
19- dir = './' ):
20+ dir = './' , logger = None ):
2021 """
2122 Build an MPAS mesh using JIGSAW with the given cell sizes as a function of
2223 latitude and longitude.
@@ -43,64 +44,70 @@ def build_spherical_mesh(cellWidth, lon, lat, earth_radius,
4344 The file name of the resulting MPAS mesh
4445
4546 plot_cellWidth : bool, optional
46- Whether to produce a plot of ``cellWidth``. If so, it will be written to
47- ``cellWidthGlobal.png``.
47+ Whether to produce a plot of ``cellWidth``. If so, it will be written
48+ to ``cellWidthGlobal.png``.
4849
4950 dir : str, optional
5051 A directory in which a temporary directory will be added with files
5152 produced during mesh conversion and then deleted upon completion.
53+
54+ logger : logging.Logger, optional
55+ A logger for the output if not stdout
5256 """
5357
54- da = xarray .DataArray (cellWidth ,
55- dims = ['lat' , 'lon' ],
56- coords = {'lat' : lat , 'lon' : lon },
57- name = 'cellWidth' )
58- cw_filename = 'cellWidthVsLatLon.nc'
59- da .to_netcdf (cw_filename )
60- if plot_cellWidth :
61- register_sci_viz_colormaps ()
62- fig = plt .figure (figsize = [16.0 , 8.0 ])
63- ax = plt .axes (projection = ccrs .PlateCarree ())
64- ax .set_global ()
65- im = ax .imshow (cellWidth , origin = 'lower' ,
66- transform = ccrs .PlateCarree (),
67- extent = [- 180 , 180 , - 90 , 90 ], cmap = '3Wbgy5' ,
68- zorder = 0 )
69- ax .add_feature (cartopy .feature .LAND , edgecolor = 'black' , zorder = 1 )
70- gl = ax .gridlines (
71- crs = ccrs .PlateCarree (),
72- draw_labels = True ,
73- linewidth = 1 ,
74- color = 'gray' ,
75- alpha = 0.5 ,
76- linestyle = '-' , zorder = 2 )
77- gl .top_labels = False
78- gl .right_labels = False
79- plt .title (
80- 'Grid cell size, km, min: {:.1f} max: {:.1f}' .format (
81- cellWidth .min (),cellWidth .max ()))
82- plt .colorbar (im , shrink = .60 )
83- fig .canvas .draw ()
84- plt .tight_layout ()
85- plt .savefig ('cellWidthGlobal.png' , bbox_inches = 'tight' )
86- plt .close ()
87-
88- print ('Step 1. Generate mesh with JIGSAW' )
89- jigsaw_driver (cellWidth , lon , lat , on_sphere = True ,
90- earth_radius = earth_radius )
91-
92- print ('Step 2. Convert triangles from jigsaw format to netcdf' )
93- jigsaw_to_netcdf (msh_filename = 'mesh-MESH.msh' ,
94- output_name = 'mesh_triangles.nc' , on_sphere = True ,
95- sphere_radius = earth_radius )
96-
97- print ('Step 3. Convert from triangles to MPAS mesh' )
98- write_netcdf (convert (xarray .open_dataset ('mesh_triangles.nc' ), dir = dir ),
99- out_filename )
58+ with LoggingContext (__name__ , logger = logger ) as logger :
59+
60+ da = xarray .DataArray (cellWidth ,
61+ dims = ['lat' , 'lon' ],
62+ coords = {'lat' : lat , 'lon' : lon },
63+ name = 'cellWidth' )
64+ cw_filename = 'cellWidthVsLatLon.nc'
65+ da .to_netcdf (cw_filename )
66+ if plot_cellWidth :
67+ register_sci_viz_colormaps ()
68+ fig = plt .figure (figsize = [16.0 , 8.0 ])
69+ ax = plt .axes (projection = ccrs .PlateCarree ())
70+ ax .set_global ()
71+ im = ax .imshow (cellWidth , origin = 'lower' ,
72+ transform = ccrs .PlateCarree (),
73+ extent = [- 180 , 180 , - 90 , 90 ], cmap = '3Wbgy5' ,
74+ zorder = 0 )
75+ ax .add_feature (cartopy .feature .LAND , edgecolor = 'black' , zorder = 1 )
76+ gl = ax .gridlines (
77+ crs = ccrs .PlateCarree (),
78+ draw_labels = True ,
79+ linewidth = 1 ,
80+ color = 'gray' ,
81+ alpha = 0.5 ,
82+ linestyle = '-' , zorder = 2 )
83+ gl .top_labels = False
84+ gl .right_labels = False
85+ plt .title (
86+ 'Grid cell size, km, min: {:.1f} max: {:.1f}' .format (
87+ cellWidth .min (),cellWidth .max ()))
88+ plt .colorbar (im , shrink = .60 )
89+ fig .canvas .draw ()
90+ plt .tight_layout ()
91+ plt .savefig ('cellWidthGlobal.png' , bbox_inches = 'tight' )
92+ plt .close ()
93+
94+ logger .info ('Step 1. Generate mesh with JIGSAW' )
95+ jigsaw_driver (cellWidth , lon , lat , on_sphere = True ,
96+ earth_radius = earth_radius , logger = logger )
97+
98+ logger .info ('Step 2. Convert triangles from jigsaw format to netcdf' )
99+ jigsaw_to_netcdf (msh_filename = 'mesh-MESH.msh' ,
100+ output_name = 'mesh_triangles.nc' , on_sphere = True ,
101+ sphere_radius = earth_radius )
102+
103+ logger .info ('Step 3. Convert from triangles to MPAS mesh' )
104+ write_netcdf (convert (xarray .open_dataset ('mesh_triangles.nc' ), dir = dir ,
105+ logger = logger ),
106+ out_filename )
100107
101108
102109def build_planar_mesh (cellWidth , x , y , geom_points , geom_edges ,
103- out_filename = 'base_mesh.nc' ):
110+ out_filename = 'base_mesh.nc' , logger = None ):
104111 """
105112 Build a planar MPAS mesh
106113
@@ -121,28 +128,31 @@ def build_planar_mesh(cellWidth, x, y, geom_points, geom_edges,
121128
122129 out_filename : str, optional
123130 The file name of the resulting MPAS mesh
131+
132+ logger : logging.Logger, optional
133+ A logger for the output if not stdout
124134 """
125- da = xarray . DataArray ( cellWidth ,
126- dims = [ 'y' , 'x' ],
127- coords = { 'y' : y , 'x' : x },
128- name = ' cellWidth' )
129- cw_filename = 'cellWidthVsXY.nc'
130- da . to_netcdf ( cw_filename )
131-
132- print ( 'Step 1. Generate mesh with JIGSAW' )
133- jigsaw_driver (
134- cellWidth ,
135- x ,
136- y ,
137- on_sphere = False ,
138- geom_points = geom_points ,
139- geom_edges = geom_edges )
140-
141- print ( 'Step 2. Convert triangles from jigsaw format to netcdf' )
142- jigsaw_to_netcdf ( msh_filename = 'mesh-MESH.msh' ,
143- output_name = 'mesh_triangles.nc' , on_sphere = False )
144-
145- print ( 'Step 3. Convert from triangles to MPAS mesh' )
146- write_netcdf ( convert ( xarray . open_dataset ( 'mesh_triangles.nc' ) ),
147- out_filename )
135+
136+ with LoggingContext ( __name__ , logger = logger ) as logger :
137+
138+ da = xarray . DataArray ( cellWidth ,
139+ dims = [ 'y' , 'x' ],
140+ coords = { 'y' : y , 'x' : x },
141+ name = 'cellWidth' )
142+ cw_filename = 'cellWidthVsXY.nc'
143+ da . to_netcdf ( cw_filename )
144+
145+ logger . info ( 'Step 1. Generate mesh with JIGSAW' )
146+ jigsaw_driver ( cellWidth , x , y , on_sphere = False ,
147+ geom_points = geom_points , geom_edges = geom_edges ,
148+ logger = logger )
149+
150+ logger . info ( 'Step 2. Convert triangles from jigsaw format to netcdf' )
151+ jigsaw_to_netcdf ( msh_filename = 'mesh-MESH.msh' ,
152+ output_name = 'mesh_triangles.nc' , on_sphere = False )
153+
154+ logger . info ( 'Step 3. Convert from triangles to MPAS mesh' )
155+ write_netcdf ( convert ( xarray . open_dataset ( 'mesh_triangles.nc' ),
156+ logger = logger ),
157+ out_filename )
148158
0 commit comments