Skip to content

Instantly share code, notes, and snippets.

@gubuntu

gubuntu/SACRS.md

Last active May 12, 2021 15:30
Show Gist options
  • Star 5 You must be signed in to star a gist
  • Fork 2 You must be signed in to fork a gist
  • Save gubuntu/6403425 to your computer and use it in GitHub Desktop.
Save gubuntu/6403425 to your computer and use it in GitHub Desktop.
Download ZIP
southern African coordinate reference systems (CRS) definitions in proj4 and WKT
Raw
SACRS.md

Use this file to get the CRS definitions you need in southern Africa, particularly South Africa, Lesotho, eSwatini, Zimbabwe and Namibia.

Most GIS software now implements the South African CRS (SACRS; incorrectly but popularly known as the "LO" system) which is based on south-oriented (or south-facing, 'west-south up') transverse mercator (TMSO). The codes for these are:

  • Hartebeesthoek94 or WGS84 datum (post 1999)
    • EPSG: all numbers from 2046 to 2055
    • ESRI (rather use EPSG): 102480 to 102489
  • Cape datum (for pre-1999 data)
    • EPSG: odd numbers from 22275 to 22293
    • ESRI (rather use EPSG): 102470 to 102479
  • The Namibian Schwarzeck CRS (on the Bessel Namibia datum) is also south-oriented: odd numbers from 29371 to 29385

If you have data with coordinates in one of these south-oriented CRSs it will (correctly) draw upside down and back to front in its native CRS but if you have on-the-fly transformation enabled, your layers will overlay your other data perfectly in whatever project CRS you choose. This feature request will result south-oriented coordinates and visualisation being more intutitive to GIS users.

You'll get south-oriented coordinates from land surveyors or CAD exports, but in practice most GIS practitioners use north-oriented (or 'north-facing' or 'east-north up') transverse mercator. So, here are some definitions you can use for north-oriented TM, as well as some for Albers, and some other CRS guidelines for southern Africa.

There are no EPSG codes for any north-oriented CRSs for South Africa or Namibia, so you have to roll your own in most software. Since Proj.7 landed in QGIS 3.x, there are some ESRI CRS definitions you can also try.

SACRS North-oriented CRS on Hartebeesthoek94 or WGS84 datum

(These are the normal TM CRS you would use in everyday GIS)

If you are working in an area up to 1 degree either side of the central meridian, this is the ideal CRS.

SACRS_NO_27 in proj4

+proj=tmerc +lat_0=0 +lon_0=27 +k=1 +x_0=0 +y_0=0 +axis=enu +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs

Note the +axis=enu, which denotes a north-oriented CRS.

SACRS_NO_27 in WKT

(what you'd get in a .prj file, or the GeoServer epsg_properties file, or recent QGIS/Proj versions)

PROJCS["SACRS_NO_27",GEOGCS["WGS 84",DATUM["unknown",SPHEROID["WGS84",6378137,298.257223563],TOWGS84[0,0,0,0,0,0,0]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",27],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1]]

To get the other TM zones, just substitute 27 with the zone you want.

These are also now available in the ESRI codes 102562 - 102568 (the (E_N) ones).

SACRS North-oriented on Cape datum

(e.g. for old South African data or Lesotho data)

TM_Cape_NO_27 in proj4

+proj=tmerc +lat_0=0 +lon_0=27 +k=1 +x_0=0 +y_0=0 +axis=enu +a=6378249.145 +b=6356514.966398753 +towgs84=-136,-108,-292,0,0,0,0 +units=m +no_defs

TM_Cape_NO_27 in WKT

             PROJCS["TM_Cape_NO_27",GEOGCS["GCS_Cape",DATUM["D_Cape",SPHEROID["Clarke_1880_Arc",6378249.145,293.4663077]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",27],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1]]

Namibian Schwarzeck example

EPSG 29375 = '+proj=tmerc +lat_0=-22 +lon_0=15 +k=1 +x_0=0 +y_0=0 +axis=wsu +ellps=bess_nam +towgs84=616,97,-251,0,0,0,0 +to_meter=1.0000135965 +no_defs'

A point defined as SRID=29375;POINT(-14242.82 86750.36) is incorrect because these are ENU ('East North Up') coordinates, not 'West South Up'.

If you want to define your geometry as 29375 you should load your raw coordinates (as they come from the surveyors) like this:

SRID=29375;POINT(14242.82 -86750.36). In other words, you do not need to modify raw surveyor coordinates. Note how the coordinates are in the same order but the signs are swapped compared to the incorrect example above.

You can fix them in place with a query like this in PostGIS so they work perfectly everywhere:

update beacons set geom = 'SRID=29375;POINT('|| -ST_x(geom) ||' '|| -ST_y(geom) ||')';

Note: If your data are in 29375 and your QGIS canvas is in 29375 everything will overlay correctly; however, it will draw upside down. To make it north-up again, just project your canvas to any other north-up CRS.

Otherwise, if you want your modified 'East North Up' coordinates to work seamlessly you need to add a custom SRID to PostGIS as follows. You can also add it explicitly to QGIS although QGIS will pick it up when it opens the layer from PostGIS (you would also need to add it to GeoServer, ArcGIS or other GIS software):

insert into spatial_ref_sys (srid,proj4text)
VALUES (40015, '+proj=tmerc +lat_0=-22 +lon_0=15 +k=1 +x_0=0 +y_0=0 +axis=enu +ellps=bess_nam +towgs84=616,97,-251,0,0,0,0 +to_meter=1.0000135965 +no_defs');

Note the only difference is +axis=enu instead of wsu.

Then you can update the SRID as follows

select updategeometrysrid('uis_beacons','geom',40015);

Now your geometries will look like this and they will work perfectly everywhere:

SRID=40015;POINT(-14242.82 86750.36)

UTM

If you are working in an area up to six degrees wide that fits into one UTM zone, then use one of the predefined UTM zones on the appropriate datum.

e.g. UTM 34S on the WGS84 datum, which is EPSG:32734

In which UTM zone does your study area fall? http://www.dmap.co.uk/utmworld.htm

Albers

If you are working in a larger area such as a whole province or country, then define your own Albers projection to suit. Here's an example appropriate for the whole of South Africa:

proj4

+proj=aea +lat_1=-24 +lat_2=-32 +lat_0=0 +lon_0=24 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs

WKT

PROJCS["Albers_SA",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Albers"],PARAMETER["standard_parallel_1",-24],PARAMETER["standard_parallel_2",-32],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",24],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1]]

Guidelines for the parameters when defining your own Albers CRS are that the central meridian should go through the centre of your area and the 1st and 2nd parallels should go one sixth down from the northern extent and one sixth up from the southern extent respectively.

Albers EPSG code for continental South Africa: 9221 (Since around 2020)

Software implementation

QGIS

Go to Settings -> Custom CRS and add a new CRS with the appropriate proj4 definition.

Try the Projestions plugin to suggest appropriate CRS for your locality.

The CRS tab in Prject Properties also graphically shows the applicable extent of most CRSs.

PostGIS

Insert a record into spatial_ref_sys with your own unique SRID and the appropriate matching proj4 and WKT (srtext) definition.

This example for an Albers CRS suitable for the whole of South Africa:

INSERT INTO spatial_ref_sys (srid,proj4text) 
VALUES (40000,'+proj=aea +lat_1=-24,241886157 +lat_2=-32,715290215 +lat_0=0 +lon_0=24,670537409 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs')

GeoServer

Add a line to <geoserver data dir>/user_projections/epsg.properties with your own unique code and a WKT definition and restart GeoServer.

These are mine for TM North facing for all the zones in SA:

40015=PROJCS["HBKNO15",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",15],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40015"]]]
40017=PROJCS["HBKNO17",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",17],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40017"]]]
40019=PROJCS["HBKNO19",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",19],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40019"]]]
40021=PROJCS["HBKNO21",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",21],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40021"]]]
40023=PROJCS["HBKNO23",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",23],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40023"]]]
40025=PROJCS["HBKNO25",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",25],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40025"]]]
40027=PROJCS["HBKNO27",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",27],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40027"]]]
40029=PROJCS["HBKNO29",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",29],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40029"]]]
40031=PROJCS["HBKNO31",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",31],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40031"]]]
40033=PROJCS["HBKNO33",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",33],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","40033"]]]

GDAL/OGR/Proj4

Either use explicit proj4 code strings for CRS definitions that are not in the EPSG database or add them to your installation's EPSG database file which in Ubuntu is /usr/share/proj/epsg. Just give them unique IDs like the ones I used above for GeoServer.

(Note, this might have changed in Proj.5, 6 and 7)

ArcGIS

Define a new projection in ArcCatalog with the appropriate WKT definition.

Note on datums

Although there are a few metres difference between the WGS84 and Hartebeethoek94 datums, for all intents and purposes in the GIS world they can be considered the same - however the EPSG does provide a datum transformation between the two: 1505 or its inverse.

If you need to transform from the Cape datum to Hartebeesthoek94, do NOT do it in your GIS. Your GIS software will apply a blanket transformation which will introduce large errors. You need instead to use specialist software and GCPs to locally transform your data because the error in the Cape datum varies hugely across South Africa.

Note added in Jan 2021: Deprecated or dropped during 2019-2020 with the releases of Proj.6 and Proj.7 and QGIS 3.x

QGIS dropped the 'ZANGI', 'SACRS_NO' CRSs, since QGIS no longer maintains its own CRS database but instead uses whatever comes with Proj, which is mainly EPSG supplemented by a few other sources like IGNF and ESRI. Furthermore, the old proj.4 CRS string definitions are now deprecated in favour of WKT descriptions. These changes have landed / will land for you depending on when QGIS was built for your OS with the versions of GDAL and PROJ that implemented them.

History: Tim Sutton and I added what I called "SACRS_NO" CRSs to QGIS for the launch of v2.0 in 2013. This is what the acronyms mean: ZANGI (ZA [South Africa] National Geospatial Information); HBK (Hartebeesthoek94); NO (north-oriented):

  • Cape datum: South African CRS : CAPE_NO_15 or ZANGI:CPNO15 (substitute any central meridian from 15 to 33)
  • Hartebeesthoek94 datum: South African CRS : HBK_NO_15 or ZANGI:HBKNO15 (substitute any central meridian from 15 to 33)

Disclaimer

Use these at your own risk.

This is a public gist - feel free to edit or comment.

@tem11010
Copy link

Thank for this. I've been trying to find this information for some time! Much appreciated. Tim

@JordaanFouche
Copy link

Thank you Tim

This is extremely useful. I tried to generate my own CRS using the WKT for 27 and got an error message like this:

image

Also I am trying to convert coordinates that are in WG27 to WGS 84 or Hartebeeshoek Lo27 and I am clearly failing. My data plots in the northern hemisphere!

I hope that you can help me.

Best regards

Jordaan

@gubuntu
Copy link
Author

gubuntu commented Jan 18, 2021

I hope your questions were answered in this thread https://lists.osgeo.org/pipermail/qgis-za-user/2021-January/000292.html

If you want to create a custom CRS for SA in QGis then you neeed to supply a WKT (preferred over proj since 3.16) definition like https://gist.github.com/gubuntu/6403425#sacrs_no_27-in-wkt. In your example you seemed to be trying to recreate 3395, which already existed.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment