@@ -4082,6 +4082,16 @@ def _safe_unit_str(unit_obj):
40824082 return str (unit_obj )
40834083
40844084
4085+ def _normalize_unit_obj (unit_obj ):
4086+ if unit_obj is None :
4087+ return None
4088+ if isinstance (unit_obj , str ):
4089+ return units (unit_obj ).units
4090+ if hasattr (unit_obj , "units" ):
4091+ return unit_obj .units
4092+ return unit_obj
4093+
4094+
40854095def _first_derivative_variable (field , delta , axis ):
40864096 """Compute first derivative with variable spacing along an axis.
40874097
@@ -4248,8 +4258,10 @@ def divergence(u, v, dx=None, dy=None, x_dim=-1, y_dim=-2, parallel_scale=None,
42484258 dx , dy = _resolve_dx_dy (v , dx = dx , dy = dy , latitude = latitude , longitude = longitude )
42494259 if dx is None or dy is None :
42504260 raise TypeError ("divergence requires dx/dy or inferable latitude/longitude coordinates" )
4251- u_arr = np .asarray (_strip (u , "m/s" ), dtype = np .float64 )
4252- v_arr = np .asarray (_strip (v , "m/s" ), dtype = np .float64 )
4261+ u_unit = _data_unit (u )
4262+ v_unit = _data_unit (v )
4263+ u_arr = np .asarray (_strip (u , str (u_unit ) if u_unit is not None else "m/s" ), dtype = np .float64 )
4264+ v_arr = np .asarray (_strip (v , str (v_unit ) if v_unit is not None else "m/s" ), dtype = np .float64 )
42534265
42544266 # Check for map projection scale factors (spherical corrections)
42554267 if parallel_scale is None and meridional_scale is None :
@@ -4266,20 +4278,26 @@ def divergence(u, v, dx=None, dy=None, x_dim=-1, y_dim=-2, parallel_scale=None,
42664278 du_dx_corr , _ , _ , dv_dy_corr = _vector_derivative_corrected (
42674279 u_arr , v_arr , dx , dy , ps , ms )
42684280 result = du_dx_corr + dv_dy_corr
4281+ if u_unit is not None or v_unit is not None :
4282+ return _wrap_derivative_like (u if u_unit is not None else v , result , units .m )
42694283 return _wrap_result_like (u , result , "1/s" )
42704284
42714285 # Variable spacing without scale factors
42724286 if _is_variable_spacing (dx ) or _is_variable_spacing (dy ) or dx_m .ndim >= 2 :
42734287 dudx = _first_derivative_variable (u_arr , dx_m , axis = - 1 )
42744288 dvdy = _first_derivative_variable (v_arr , dy_m , axis = - 2 )
42754289 result = dudx + dvdy
4290+ if u_unit is not None or v_unit is not None :
4291+ return _wrap_derivative_like (u if u_unit is not None else v , result , units .m )
42764292 return _wrap_result_like (u , result , "1/s" )
42774293
42784294 # Uniform grid: fast Rust path
42794295 dx_val = float (dx_m .mean ()) if dx_m .ndim > 0 else float (dx_m )
42804296 dy_val = float (dy_m .mean ()) if dy_m .ndim > 0 else float (dy_m )
42814297 if u_arr .ndim == 2 :
42824298 result = np .asarray (_calc .divergence (np .ascontiguousarray (u_arr ), np .ascontiguousarray (v_arr ), dx_val , dy_val ))
4299+ if u_unit is not None or v_unit is not None :
4300+ return _wrap_derivative_like (u if u_unit is not None else v , result , units .m )
42834301 return _wrap_result_like (u , result , "1/s" )
42844302 lead_shape = u_arr .shape [:- 2 ]
42854303 ny , nx = u_arr .shape [- 2 :]
@@ -4290,6 +4308,8 @@ def divergence(u, v, dx=None, dy=None, x_dim=-1, y_dim=-2, parallel_scale=None,
42904308 result .reshape ((- 1 , ny , nx ))[idx ] = np .asarray (
42914309 _calc .divergence (np .ascontiguousarray (u_flat [idx ]), np .ascontiguousarray (v_flat [idx ]), dx_val , dy_val )
42924310 ).reshape ((ny , nx ))
4311+ if u_unit is not None or v_unit is not None :
4312+ return _wrap_derivative_like (u if u_unit is not None else v , result , units .m )
42934313 return _wrap_result_like (u , result , "1/s" )
42944314
42954315
@@ -4655,22 +4675,28 @@ def _cross_section_distance_coords(cross):
46554675 lon ,
46564676 lat ,
46574677 )
4678+ x_vals = np .asarray (distance * np .sin (np .deg2rad (forward_az )), dtype = np .float64 )
4679+ y_vals = np .asarray (distance * np .cos (np .deg2rad (forward_az )), dtype = np .float64 )
46584680 x = xr .DataArray (
4659- units . Quantity ( distance * np . sin ( np . deg2rad ( forward_az )), "meter" ) ,
4681+ x_vals ,
46604682 coords = lon_coord .coords ,
46614683 dims = lon_coord .dims ,
4684+ attrs = {"units" : "meter" },
46624685 )
46634686 y = xr .DataArray (
4664- units . Quantity ( distance * np . cos ( np . deg2rad ( forward_az )), "meter" ) ,
4687+ y_vals ,
46654688 coords = lat_coord .coords ,
46664689 dims = lat_coord .dims ,
4690+ attrs = {"units" : "meter" },
46674691 )
46684692 return x , y
46694693 if "x" in coords and "y" in coords :
46704694 x_coord = coords ["x" ]
46714695 y_coord = coords ["y" ]
4672- x = xr .DataArray (_dataarray_unit_array (x_coord , "meter" ), coords = x_coord .coords , dims = x_coord .dims )
4673- y = xr .DataArray (_dataarray_unit_array (y_coord , "meter" ), coords = y_coord .coords , dims = y_coord .dims )
4696+ x_vals = np .asarray (_dataarray_unit_array (x_coord , "meter" ).to ("meter" ).magnitude , dtype = np .float64 )
4697+ y_vals = np .asarray (_dataarray_unit_array (y_coord , "meter" ).to ("meter" ).magnitude , dtype = np .float64 )
4698+ x = xr .DataArray (x_vals , coords = x_coord .coords , dims = x_coord .dims , attrs = {"units" : "meter" })
4699+ y = xr .DataArray (y_vals , coords = y_coord .coords , dims = y_coord .dims , attrs = {"units" : "meter" })
46744700 return x , y
46754701 raise AttributeError ("Sufficient horizontal coordinates not defined." )
46764702
@@ -4680,16 +4706,18 @@ def _cross_section_latitude_coord(cross):
46804706 if "latitude" in coords :
46814707 latitude = coords ["latitude" ]
46824708 return xr .DataArray (
4683- _dataarray_unit_array (latitude , "degree" ),
4709+ np . asarray ( _dataarray_unit_array (latitude , "degree" ). to ( "degree" ). magnitude , dtype = np . float64 ),
46844710 coords = latitude .coords ,
46854711 dims = latitude .dims ,
4712+ attrs = {** dict (latitude .attrs ), "units" : latitude .attrs .get ("units" , "degree" )},
46864713 )
46874714 if "lat" in coords :
46884715 latitude = coords ["lat" ]
46894716 return xr .DataArray (
4690- _dataarray_unit_array (latitude , "degree" ),
4717+ np . asarray ( _dataarray_unit_array (latitude , "degree" ). to ( "degree" ). magnitude , dtype = np . float64 ),
46914718 coords = latitude .coords ,
46924719 dims = latitude .dims ,
4720+ attrs = {** dict (latitude .attrs ), "units" : latitude .attrs .get ("units" , "degree" )},
46934721 )
46944722 if hasattr (cross , "metpy" ):
46954723 try :
@@ -4701,8 +4729,12 @@ def _cross_section_latitude_coord(cross):
47014729 inverse = True ,
47024730 radians = False ,
47034731 )[1 ]
4704- return xr .DataArray (units .Quantity (latitude , "degrees_north" ), coords = y .coords ,
4705- dims = y .dims )
4732+ return xr .DataArray (
4733+ np .asarray (latitude , dtype = np .float64 ),
4734+ coords = y .coords ,
4735+ dims = y .dims ,
4736+ attrs = {"units" : "degrees_north" },
4737+ )
47064738 except Exception :
47074739 pass
47084740 raise AttributeError ("Latitude coordinates are required for absolute_momentum." )
@@ -5114,7 +5146,12 @@ def absolute_momentum(u, v=None, index="index"):
51145146 f = coriolis_parameter (latitude )
51155147 x , y = _cross_section_distance_coords (norm_wind )
51165148 distance_q = np .hypot (_dataarray_unit_array (x , "meter" ), _dataarray_unit_array (y , "meter" ))
5117- distance = xr .DataArray (distance_q , coords = x .coords , dims = x .dims )
5149+ distance = xr .DataArray (
5150+ np .asarray (distance_q .to ("meter" ).magnitude , dtype = np .float64 ),
5151+ coords = x .coords ,
5152+ dims = x .dims ,
5153+ attrs = {"units" : "meter" },
5154+ ).metpy .quantify ()
51185155 _ , distance = xr .broadcast (norm_wind , distance )
51195156 result = norm_wind + f * distance
51205157 return result .metpy .convert_units ("m/s" ) if hasattr (result , "metpy" ) else result
@@ -5142,12 +5179,14 @@ def coriolis_parameter(latitude):
51425179 lat_arr = np .asarray (lat_stripped , dtype = np .float64 )
51435180 result = 2.0 * 7.292115e-5 * np .sin (np .deg2rad (lat_arr ))
51445181 if _is_dataarray_like (latitude ):
5182+ attrs = dict (latitude .attrs )
5183+ attrs ["units" ] = "1/s"
51455184 return xr .DataArray (
5146- np .asarray (result , dtype = np .float64 ) * units ( "1/s" ) ,
5185+ np .asarray (result , dtype = np .float64 ),
51475186 coords = latitude .coords ,
51485187 dims = latitude .dims ,
5149- attrs = dict ( latitude . attrs ) ,
5150- )
5188+ attrs = attrs ,
5189+ ). metpy . quantify ()
51515190 return _wrap_result_like (latitude , result , "1/s" )
51525191
51535192
@@ -5216,7 +5255,8 @@ def curvature_vorticity(u, v, dx=None, dy=None, x_dim=-1, y_dim=-2, parallel_sca
52165255 longitude = longitude ,
52175256 crs = crs ,
52185257 )
5219- return (u * u * dvdx - v * v * dudy - v * u * dudx + u * v * dvdy ) / (u ** 2 + v ** 2 )
5258+ with np .errstate (divide = "ignore" , invalid = "ignore" ):
5259+ return (u * u * dvdx - v * v * dudy - v * u * dudx + u * v * dvdy ) / (u ** 2 + v ** 2 )
52205260
52215261
52225262def inertial_advective_wind (u , v , u_geostrophic , v_geostrophic , dx = None , dy = None ,
@@ -5362,7 +5402,8 @@ def shear_vorticity(u, v, dx=None, dy=None, x_dim=-1, y_dim=-2, parallel_scale=N
53625402 longitude = longitude ,
53635403 crs = crs ,
53645404 )
5365- return (v * u * dudx + v * v * dvdx - u * u * dudy - u * v * dvdy ) / (u ** 2 + v ** 2 )
5405+ with np .errstate (divide = "ignore" , invalid = "ignore" ):
5406+ return (v * u * dudx + v * v * dvdx - u * u * dudy - u * v * dvdy ) / (u ** 2 + v ** 2 )
53665407
53675408
53685409def shearing_deformation (u , v , dx = None , dy = None , x_dim = - 1 , y_dim = - 2 ,
@@ -6449,16 +6490,16 @@ def _gradient_spacing(position, axis_size, mode):
64496490
64506491def _data_unit (data ):
64516492 if hasattr (data , "units" ):
6452- return data .units
6493+ return _normalize_unit_obj ( data .units )
64536494 try :
64546495 unit_attr = data .metpy .units
64556496 if str (unit_attr ) != "dimensionless" :
6456- return unit_attr
6497+ return _normalize_unit_obj ( unit_attr )
64576498 except Exception :
64586499 pass
64596500 attrs = getattr (data , "attrs" , {})
64606501 if isinstance (attrs , dict ) and "units" in attrs :
6461- return units (attrs ["units" ])
6502+ return _normalize_unit_obj (attrs ["units" ])
64626503 return None
64636504
64646505
0 commit comments