search

csyhuang / hn2016_falwa

Python library for computing Finite-Amplitude Local Wave Activity Diagnostics from climate data.
https://hn2016-falwa.readthedocs.io/en/latest/
MIT License
32 stars 16 forks source link

A potential bug in the 'u_baro' calculation of FALWA v.1.2.0 #108

Closed sandrolubis closed 8 months ago

sandrolubis commented 8 months ago

Hi all, after updating to FALWA v1.2.0, I spotted a potential bug in u_baro calculation where the output mirrors - i.e., SH = NH. I checked the 'interpolated_u' data; they are correct. Only 'u_baro' shows this issue. See the attached plot.

ubaro_falwa120

chpolste commented 8 months ago

I can reproduce the problem with QGFieldNHN22. QGFieldNH18 does not seem to be affected.

I think the loop in

https://github.com/csyhuang/hn2016_falwa/blob/e47f685846a0186682ba6d473174d6ab8e4127be/falwa/f90_modules/compute_flux_dirinv.f90#L190-L193

needs to be adapted to the value of is_nhem as in the loop earlier in the function. Maybe remove

https://github.com/csyhuang/hn2016_falwa/blob/e47f685846a0186682ba6d473174d6ab8e4127be/falwa/f90_modules/compute_flux_dirinv.f90#L181-L188

and replace everything with

    if (is_nhem) then
      do j = jstart,jend
        ubaro(:,j) = ubaro(:,j)+uu(:,j+nd-1,k)*exp(-zk/h)*dc
        urefbaro(j) = urefbaro(j)+uref(j-jb,k)*exp(-zk/h)*dc
      enddo
    else
      do j = jstart,jend
        ubaro(:,j) = ubaro(:,j)+uu(:,j,k)*exp(-zk/h)*dc
        urefbaro(j) = urefbaro(j)+uref(j,k)*exp(-zk/h)*dc
      enddo
    endif

@csyhuang, do you think this works with the indices? It seems to resolve the issue in my test notebook, but I'm worried about an off-by-one error...

csyhuang commented 8 months ago

@sandrolubis Thanks for reporting the bug! @chpolste Thanks for figuring out the cause of this issue! I think your solution works (maybe also define explicitly jstart and jend in the loop?)! Would you like to submit a bugfix PR (and update version number to v1.2.1)? If not, I can do it this weekend after work. Thanks again 🙏 😊