wrf-python icon indicating copy to clipboard operation
wrf-python copied to clipboard

Issue with vertcross function when using Fortran-order data from .npy files

Open Wang5230 opened this issue 1 year ago • 2 comments

I encountered an issue when using the vertcross function from the wrf-python package to perform vertical interpolation on data stored in a .npy file. Specifically, when I load the .npy file , the data is stored in Fortran-order (column-major) format by default.

The vertcross function seems to require C-order (row-major) data storage format to work correctly. When passing Fortran-order data, I received a ValueError: invalid shape for argument 1 (got (38,), expected (360,)), which suggests a misinterpretation of the data shape due to the storage format mismatch.

Steps to Reproduce:

  1. Load vorticity data counted by wrfout (e.g., vorticity) from a .npy file:

    vorticity = np.load("vorticity.npy")
    
  2. Attempt to use vertcross to perform vertical interpolation:

    vs_cross = vertcross(vorticity, z, ...)
    
  3. Encounter the error: ValueError: invalid shape for argument 1 (got (38,), expected (360,)) The vertcross function should work correctly when the data is in Fortran-order (which is the default when reading .npy files). If the function requires C-order, it would be helpful to include a clear message in the documentation or provide an automatic conversion inside the function.

solution:

Manually converting the data to C-order format using np.ascontiguousarray() resolves the issue:

vorticity = np.ascontiguousarray(vorticity)
vs_cross = vertcross(vorticity, z, ...)

It would be helpful if the vertcross function could automatically handle Fortran-order data or provide better error handling or documentation related to data format expectations.

Wang5230 avatar Mar 11 '25 08:03 Wang5230

Behavior noted here: https://github.com/NCAR/wrf-python/blob/eec6637df55b05bc15112d1fb3c5c3917f26872a/src/wrf/extension.py#L178

Switching both arrays to fortran-contiguous z-y-x produces an array of nans:

>>> from wrf import vertcross, CoordPair
>>> vorticity = np.zeros((61, 184, 251), dtype="f4", order="f")
>>> vertical = vorticity.copy()
>>> vertical[:, :, :] = np.arange(61)[:, np.newaxis, np.newaxis]
>>> center = CoordPair(x=120, y=90)
>>> vertcross(np.asfortranarray(vertical), np.asfortranarray(vertical), pivot_point=center, angle=45)
<xarray.DataArray 'field3d_cross' (dim_0: 100, dim_1: 259)> Size: 104kB
array([[nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       ...,
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan]], dtype=float32)
Dimensions without coordinates: dim_0, dim_1
Attributes:
    orientation:    (30.0, 0.0) to (213.0, 183.0) ; center=CoordPair(x=120, y...
    missing_value:  9.969209968386869e+36
    _FillValue:     9.969209968386869e+36

Do you want to make a PR promoting that comment to a note in the docstring and propagate that note up the call chain? The alternative would be to ensure the array is C-contiguous at the relevant point, but I'm not sure where that is.

DWesl avatar Mar 11 '25 22:03 DWesl

Thanks for your response! I'll take some time to look deeper into the code to determine whether it's better to submit a PR that adds a note in the docstring or to modify the code to ensure the array is C-contiguous automatically. This might take some time, but I'll follow up once I have a better understanding.

Wang5230 avatar Mar 12 '25 13:03 Wang5230