TEOS-10 / GibbsSeaWater.jl

Gibbs-SeaWater (GSW) Oceanographic Toolbox in Julia
http://www.teos-10.org
Other
23 stars 3 forks source link

Manual vs automatic modile update #4

Open ax1ine opened 4 years ago

ax1ine commented 4 years ago

I found some differences between the latest C header and Julia wrappers. For example, these functions:

extern double gsw_z_from_p(double p, double lat, double geo_strf_dyn_height, double sea_surface_geopotential);

extern double gsw_p_from_z(double z, double lat, double geo_strf_dyn_height, double sea_surface_geopotential); go as

function gsw_z_from_p(p::Cdouble, lat::Cdouble) ccall((:gsw_z_from_p, libgswteos), Cdouble, (Cdouble, Cdouble), p, lat) end

function gsw_p_from_z(z::Cdouble, lat::Cdouble) ccall((:gsw_p_from_z, libgswteos), Cdouble, (Cdouble, Cdouble), z, lat) end

Maybe there are some more inconsistent functions. The questions is: how are we going to update the module? Are we going to do it automatically via CLang.jl, or manually?

Would be nice to add comments to the code (to make the module compliant to Julia documentations standards) but I'm not sure how to be in the case of automatic updates.

Alexander-Barth commented 4 years ago

I just checked two package from JuliaGeo (LibSpatialIndex.jl and Proj4.jl) and it seems that they also start with an automatically generate wrapper file but then, as far as I can tell, update the wrapper file manually. My impression is that changes in the C API should be rather rare and adding documentation (as you mentioned), are much easier if we do not regenerate the files automatically. Also there are probably some changes to the Julia API necessary (e.g. adding a bang ! at the end of function who modify their arguments).

ax1ine commented 4 years ago

yeah, the C API is rather stable and updated one in a few months (according to the repo activity). So if there are no objections I'll add documentation to the module (and check if it's consistent with the latest C API).

dankelley commented 4 years ago

Here are a few notes from GSW-R, in case they are of use:

Alexander-Barth commented 4 years ago

I prefer SA and TC too for absolute salinity and conservative temperature. It is a good idea to use the same capitalization than the R-API (also used in matlab and Python, as far as I can tell).

The @templates in R are really nice. Unfortunately, there is no direct correspondence in Julia as far as I know.

PaulMBarker commented 4 years ago

Please keep the standard abbreviations of SA and CT for Absolute Salinity and Conservative Temperature.

hafez-ahmad commented 4 years ago

I am using Julia 1.4.2

C=[45.8;34.7] T=[28.9;22.8] P=[10.0;50.0] sp=gsw_sp_from_c(C,T,P) ERROR: MethodError: Cannot convert an object of type Array{Float64,1} to an object of type Float64 Closest candidates are: convert(::Type{T}, ::T) where T<:Number at number.jl:6 convert(::Type{T}, ::Number) where T<:Number at number.jl:7 convert(::Type{T}, ::Base.TwicePrecision) where T<:Number at twiceprecision.jl:250 ... Stacktrace: [1] cconvert(::Type{T} where T, ::Array{Float64,1}) at .\essentials.jl:390 [2] gsw_sp_from_c(::Array{Float64,1}, ::Array{Float64,1}, ::Array{Float64,1}) at C:\Users\hafez.julia\packages\GibbsSeaWater\NcNJk\src\gen_gswteos10.jl:622 [3] top-level scope at none:0

dankelley commented 4 years ago

I don't know about this problem with the Julia version, but in case it helps (as a test, when you get Julia working) the following shows an R test with your numbers.

> library(gsw)
> C<-c(45.8,34.7)
> T<-c(28.9,22.8)
> P<-c(10.0,50.0)
> gsw_SP_from_C(C,T,P)
[1] 27.28890 22.88687
hafez-ahmad commented 4 years ago

I am using Julia 1.4.2, I tried in python and R, in that cast gsw worked fine. but Julia gives these errors.

may be I am trying it in the wrong way

kouketsu commented 4 years ago

@hafez-ahmad Please check the essential use of broadcasting in Julia. So for array arguments, the following script may work fine (please notice "." after the function name).

sp=gsw_sp_from_c.(C,T,P)

In the future, we may add alias function for broadcasting :).

Alexander-Barth commented 4 years ago

I was about to answer the same thing :-) I have added this example the README for our future new users: https://github.com/TEOS-10/GibbsSeaWater.jl/commit/cf9e53ec4f8ec64e9f1aa357ae181255eacaaf5b

kouketsu commented 4 years ago

Thanks, @Alexander-Barth.

hafez-ahmad commented 4 years ago

thank you so much @kouketsu

hafez-ahmad commented 4 years ago

sorry again, I am getting error , [julia1.4.2]

gsw_geo_strf_dyn_height

ct=[28.2,22.05] # conservative temperture sa=[28.8,22.7] # absolute salinity p=[10.0,50.0] # sea pressure p_ref=500.0 gsw_geo_strf_dyn_height.(sa,ct,p,p_ref)

ERROR: MethodError: no method matching gsw_geo_strf_dyn_height(::Float64, ::Float64, ::Float64, ::Float64) Closest candidates are: gsw_geo_strf_dyn_height(::Any, ::Any, ::Any, ::Any, ::Any, ::Any) at C:\Users\hafez.julia\packages\GibbsSeaWater\NcNJk\src\gen_gswteos10.jl:234 Stacktrace: [1] _broadcast_getindex_evalf at .\broadcast.jl:631 [inlined] [2] _broadcast_getindex at .\broadcast.jl:604 [inlined] [3] getindex at .\broadcast.jl:564 [inlined] [4] copy at .\broadcast.jl:854 [inlined] [5] materialize(::Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Nothing,typeof(gsw_geo_strf_dyn_height),Tuple{Array{Float64,1},Array{Float64,1},Array{Float64,1},Float64}}) at .\broadcast.jl:820 [6] top-level scope at none:0

kouketsu commented 4 years ago

@hafez-ahmad I have added wrapper functions for dynamic height calculation. Please get the master repository (pkg> add GibbsSeaWater#master), and try like the following script.

ct=[28.2,22.05] # conservative temperture
sa=[28.8,22.7] # absolute salinity
p=[10.0,50.0] # sea pressure
p_ref=500.0
gsw_geo_strf_dyn_height(sa,ct,p,p_ref)

Note that the definition of the function is "gsw_ge_strf_dyn_height(sa::AbstractArray, ct::AbstractArray, p::AbstractArray, p_ref::Real)", so please use it without dot (".").

hafez-ahmad commented 4 years ago

I am still error ct=[28.2,22.05] # conservative temperture sa=[28.8,22.7] # absolute salinity p=[10.0,50.0] # sea pressure p_ref=500.0 gsw_geo_strf_dyn_height(sa,ct,p,p_ref)

ERROR: MethodError: no method matching gsw_geo_strf_dyn_height(::Array{Float64,1}, ::Array{Float64,1}, ::Array{Float64,1}, ::Int64) Closest candidates are: gsw_geo_strf_dyn_height(::Any, ::Any, ::Any, ::Any, ::Any, ::Any) at C:\Users\hafez.julia\packages\GibbsSeaWater\vq3GE\src\gen_gswteos10.jl:234 Stacktrace: [1] top-level scope at none:0

kouketsu commented 4 years ago

Maybe, there are two issues.

  1. the way to introduce the package. In this case, before installing the package, we should remove the old one.

    pkg> rm GibbsSeaWater
    pkg> gc
    pkg> add GibbsSeaWater#master
  2. My example is not appropriate. The original subroutine of "gsw_geo_strf_dyn_height" is buggy (please see GSW-C source) mainly due to interpolation. They recommend the alternative of "gsw_geo_strf_dyn_height_1". As the results with a little bit longer profile seemed to be successful even with "gsw_geo_strf_dyn_height", I should have raised the example with the longer profile / "gsw_geo_strf_dyn_height_1". Please try this.

    SA = [34.7118; 34.8915; 35.0256; 34.8472; 34.7366; 34.7324;]
    CT = [28.8099; 28.4392; 22.7862; 10.2262;  6.8272;  4.3236;]
    p =  [     10.0;      50;     125;     250;     600;    1000;]
    p_ref = 1000.0
    dh = GibbsSeaWater.gsw_geo_strf_dyn_height(SA, CT, p, p_ref)