JuliaClimate / ClimateTools.jl

Climate science package for Julia
https://juliaclimate.github.io/ClimateTools.jl/stable/
Other
116 stars 18 forks source link

Spatial averaging (in plot function) #109

Open ikroener opened 4 years ago

ikroener commented 4 years ago

Hi all,

I am not sure if this is the right repository since I just found the idea of outsourcing the topic of mapping/plotting to an own package (ClimateMapping).

Nevertheless, if used for Gaussian lon x lat grids the plot function needs to account for converging grid point at higher latitudes when doing a spatial average. When used to plot a global mean, the current function overemphasises polar regions. I already have a latitudinal weighted version, but wanted to open the issue before making a pull request. Furthermore, I hesitate since my version is more or less only valid on Gaussian grids.

Any ideas, any thoughts on the issue of spatial averaging in the plot function. The area-weighted averaging could also be implemented as a separate function, but then using the plot function for a global field still gives an invalid figure.

Best,

Balinus commented 4 years ago

Hello!

Thanks for this issue. The spatial averaging is indeed not a very good option (it was mostly for diagnostic tool). I also usually work on small spatial domain.

I'd be very happy if you could push a PR. From there, we could think about generalizing the approach. Perhaps looking at GeoStats might be a good idea? @juliohm

Thanks again!

juliohm commented 4 years ago

Thank you @Balinus for pinging and @ikroener for raising this spherical coordinates averaging issue. If you have experience with spherical coordinates and would like to contribute methods, please feel free to get in touch. There is a lot to be improved in GeoStats.jl but this use case is certainly included. One thing that is currently possible is to use a distance like Haversine in the estimation/simulation solvers.

Regarding area averages, there is a more general issue on my todo list: https://github.com/juliohm/GeoStats.jl/issues/40

ikroener commented 4 years ago

I finally started to setup a flux conservative remapping scheme FluxConservative, which is in its extreme case an area weighted global mean - for the moment it is the only thing it can do as well as it is limited to regular (lon x lat) grids. But I have great plans for it :) The - to be - package aims to implement the scheme as a solver for the GeoStats.EstimationProblem, so we could extend the ClimateTools.regrid function with area-weighted remapping. As far as I can see at the moment is only nearest neighbours interpolation on inverse distance possible?

@juliohm In your general todo issue you only mention that you need a fitting data structure. Does that mean you already have a suitable (julia-)method to estimate areas on a sphere? Spontaneously only something like Voronoi-Diagrams come into my mind. There are some packages in other languages, but I am not sure about the current status in julia and if they are out of the box suitable for spherical coordinates.

juliohm commented 4 years ago

And EstimationSolver would be great @ikroener , please let me know if you have any questions regarding the API and current dev tools.

We don't have the methods yet that take the volumes into account. Kriging-based solvers can deal with volumes well, but we first need a brainstorm about the data structures. After we have efficient data structures for representing the volumes (e.g. areas in 2D) we can devise efficient methods. It is great to have you in the loop, spherical coordinates is something that I would need to investigate more in order to have a good design.

In general, the data and domain types in GeoStats.jl will have two different use cases. In the first use case, they will be used as "zero-volume" collections of points with a specific spatial structure. For example, a regular grid is a collection of points with regular spacing. Most solvers, including the LearningSolvers integrated with MLJ.jl deal with this use case. The other use case takes volumes into account. So we need a trait system where we can query the volume of a cell in the regular grid volume(grid, i). In this case, spatial algorithms and Kriging-based solvers could call these functions to have better estimates.

Balinus commented 4 years ago

The - to be - package aims to implement the scheme as a solver for the GeoStats.EstimationProblem, so we could extend the ClimateTools.regrid function with area-weighted remapping. As far as I can see at the moment is only nearest neighbours interpolation on inverse distance possible?

Indeed. I transitioned the regrid method from using Python Scipy to GeoStats some months ago. The ultimate aim was/is to fully incorporate GeoStats solvers into the regrid function. As you can see though, InvDist is hardcoded :/ It served the purpose I had for the time being, but there is a lot of benefits for the climate community if we can call the GeoStats methods from the regrid.

Perhaps the most transparent to the end-user would be to simply pass solver as a regrid argument.

solver = Kriging(
  :precipitation => (variogram=GaussianVariogram(range=35.),)
)
R = ClimateTools.regrid(A, B, solver=solver)

I'm working right now on some bias correction methods for precipitation extremes (bias_ext branch). I'll try to do a version of regrid that is more parametric after.