Support integration of spherical and globe meshes

[finnlindgren]

From https://groups.google.com/g/r-inla-discussion-group/c/1f6Lj6ZbbEc/m/AtpuVhGxAgAJ :

Our samplers is a SpatialPolygonsDataFrame (called “sp”) in lat-lon coordinates, and the domain is an inla.mesh object (called mesh1) on the sphere (S^2) in Cartesian geocentric xyz coordinates. The reason that our mesh is in Cartesian xyz is we are modelling global processes, so our mesh has to be on S^2. We find we have to supply loc= Cartesian, geodetic xyz coordinates to inla.mesh.2d() to make meshes with full global coverage (i.e. to have longitude-wrapping and equal area). The inla.mesh.2d function seems to project lon-lat locations onto R^2.

When I ran ipoints(samplers=sp, domain=mesh1, group=”index”), I got the error

sp::over(samplers, integ_sp, returnList = TRUE) : identicalCRS(x, y) is not TRUE

Then I tried to use inla.spTransform(mesh1, inla.CRS("longlat")) to make mesh1 have the same CRS as the SpatialPolygons, but I got the error 'crs0' is an invalid coordinate reference object.

Part of reply:

The ipoints function doesn't have full support for global meshes yet,
as it relies on conversion to a flat equal-area projection to compute
the integration weights.
I hope to soon have a new ipoints version to properly support this use case.

It's possible this issue can't be fully solved until we get more of
the more general sf support working.

A potential workaround would be to make an integration scheme for the
entire mesh, with the mesh nodes as integration points; something like
this:

ips_sphere <- SpatialPointsDataFrame(mesh$loc, data=data.frame(
    weight = inla.mesh.fem(mesh)$va * earth.radius^2
  ), crs = fm_CRS("sphere"))
ips_globe <- fm_spTransform(ips_sphere, fm_CRS("globe"))

and then filter out the points outside the polygon, with a direct call
to sp::over. That way, the mesh wouldn't need to be involved in the
coordinate transformations.

For full integration accuracy along the sampler boundary, this step
needs to be done before aggregating the integration points, but it
might be sufficient as a temporary workaround.

The full version of the "workaround" is also the likely full solution, so something might be doable in the short term as well.