[R-sig-Geo] How to define a local equidistant CRS

[Micha Silver]

The RAINLINK package [1] derives rain rate from the attenuation of
signal strength between microwave towers. Data typically comes from
cellular network providers, and the location of the towers is
typically given in longitude/latitude.
Among the functions in RAINLINK, some require to get the attenuation
for each microwave link (pair of towers), from those towers that are
nearby, within, say 15 km. This means that the tower locations need to
be transformed to an equidistant Cartesian CRS.

In the current version of RAINLINK, this was done as follows:
An sp object was prepared from the tower locations, with CRS
"EPSG:4326". The average latitude and longitude values were obtained
as "YMiddle" and "XMiddle". Then a new, one-off azimuthal equidistant
CRS was defined as:
# Set projection string
projstring <- paste("+proj=aeqd +a=6378.137 +b=6356.752 +R_A +lat_0=",
    YMiddle, " +lon_0=",
    XMiddle," +x_0=0 +y_0=0",sep="")

This worked until the recent advances in proj and gdal. Now with proj
>=6.x and gdal >=3.x the above definition causes an error:

Error in .spTransform_Polygon(input[[i]], to_args = to_args, from_args
= from_args,  :
  coordinate operation creation from WKT failed: generic error of unknown origin
In addition: Warning messages:
1: In showSRID(uprojargs, format = "PROJ", multiline = "NO",
prefer_proj = prefer_proj) :
  Discarded ellps unknown in CRS definition: +proj=aeqd +a=6378.137
+b=6356.752 +R_A +lat_0=30.4 +lon_0=35.1 +x_0=0 +y_0=0 +type=crs
2: In showSRID(uprojargs, format = "PROJ", multiline = "NO",
prefer_proj = prefer_proj) :
  Discarded datum unknown in CRS definition

If I change the CRS definition as follows:
projstring <- paste("+proj=aeqd +R_A +lat_0=", YMiddle,
                    " +lon_0=", XMiddle,
                      " +x_0=0 +y_0=0 +ellps=WGS84",sep="")

Then using proj 6.3.2 we get a warning:
"In showSRID(uprojargs, format = "PROJ", multiline = "NO", ... :
Discarded datum Unknown based on WGS84 ellipsoid in CRS definition"

whereas on a newer setup, with "Loaded PROJ runtime: Rel. 7.2.1" the
above definition is accepted silently.

Our question is what is the reliable way to define a local Cartesian
CRS which would work anywhere in the world, and preferably be
equidistant. Should we just find the local UTM zone and use that?  UTM
is not equidistant, but it is ubiquitous, and if the cellular network
is not too broad, I assume that distances will be fairly accurate.

Best regards,
Micha

[1] https://github.com/overeem11/RAINLINK

-- 
Micha Silver
Ben Gurion Univ
Sde-Boker Remote Sensing Lab
cell: +972 (52) 3665918