Re: [OSGeo-Discuss] Parsing well-known text
On Sunday, 6 September 2020 18:31:31 EDT Bruce Bannerman wrote: > Hello Pierre, > > It is dangerous to assume what the coordinates in a data set represent. It > may not have implications for your particular use, but could have serious > consequences for other uses downstream. What would be other uses downstream? > A good approach is to consult the discovery metadata that describes the > dataset that you are using. It should help you to interpret the data. Any > good data provider will also provide discovery metadata, typically in ISO > 19115 format (or derivative). What's the relationship between WKT and ISO 19115? The point clouds I have have a description of the coordinate system in WKT or no description of the coordinate system. > If such metadata is not available, ask the data provider to better define > their data, in your case, the spatial reference system used. If they can’t, > or won’t provide this information then I would question the need to use > that data set and look for an alternative, as the data quality is suspect. The data provider may be the surveyor himself, who has flown the site with a drone. > It is important that blind assumptions of what an unqualified coordinate in > a spatial reference system means is not built into software that will be > used by others without explicitly stating that this is what is being done. The formats PerfectTIN exports to do not allow including WKT, except maybe LandXML (which I haven't looked for WKT in, but I know that one can specify the units in) and Carlson TIN (it's proprietary and undocumented, and I don't know what's in the header). The PerfectTIN file format does not have a place to put WKT, but I can add it. What I'm planning to do in Wolkenbase is this: If there's a WKT record, use it for units, and copy it to the output. If the units in WKT differ from those in Wolkenbase, convert the measurements to the unit set in Wolkenbase and edit the WKT. If there is no WKT, assume the units are as set in Wolkenbase, and add a minimal WKT which specifies only the units. If the point cloud consists of several files (which is quite common), and the variable-length records in different files disagree, I'm not sure what to do. PerfectTIN currently reads variable-length records but ignores them. I'm planning to make it use the units in WKT, if any, instead of the unit set in PerfectTIN, when converting coordinates when reading a file. Pierre -- The Black Garden on the Mountain is not on the Black Mountain. ___ Discuss mailing list Discuss@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/discuss
Re: [OSGeo-Discuss] Parsing well-known text
Hello Pierre, It is dangerous to assume what the coordinates in a data set represent. It may not have implications for your particular use, but could have serious consequences for other uses downstream. A good approach is to consult the discovery metadata that describes the dataset that you are using. It should help you to interpret the data. Any good data provider will also provide discovery metadata, typically in ISO 19115 format (or derivative). If such metadata is not available, ask the data provider to better define their data, in your case, the spatial reference system used. If they can’t, or won’t provide this information then I would question the need to use that data set and look for an alternative, as the data quality is suspect. It is important that blind assumptions of what an unqualified coordinate in a spatial reference system means is not built into software that will be used by others without explicitly stating that this is what is being done. Kind regards, Bruce > On 7 Sep 2020, at 02:53, Pierre Abbat wrote: > > On Wednesday, 2 September 2020 05:41:40 EDT Even Rouault wrote: >> As far as I know, you can't create a CRS WKT with just unit information. The >> most minimal content that validates the WKT1 grammar would be something >> like: >> >> LOCAL_CS["unspecified CRS", >>LOCAL_DATUM["unspecified datum",2000], >>UNIT["metre",1, >>AUTHORITY["EPSG","9001"]], >>AXIS["Easting",EAST], >>AXIS["Northing",NORTH]] > > That doesn't seem to have an elevation unit; should I assume that they're > both > meters? > >> Yes, PROJ >= 6 has support for parsing and creating WKT in several versions >> of the WKT standard. >> >> See proj_create_from_wkt() at >> https://proj.org/development/reference/functions.html, and all other >> proj_ getters. >> >> For creation of WKT, you might need the more advanced functions of >> https://github.com/OSGeo/PROJ/blob/master/src/proj_experimental.h , before >> exporting with proj_as_wkt() > > That's overkill for what I'm doing. I don't need to make a projection, all I > need to know is the units. I don't need to create a WKT except in the case > that the original cloud doesn't have one and I'm adding one to specify the > units. > > The code to parse the WKT is in YACC, which I don't know how to work with, > and > I don't need all the tokens. Treating them as opaque atoms and looking for > "UNIT" should work. Where can I find the grammar of the language of which WKT > is an instance? > > Adding any library multiplies the complexity of building and packaging on > Windows. PROJ is 79 times as big as Wolkenbase so far and 17 times as big as > PerfectTIN. Even Bezitopo, which handles projections, but only conformal > ones, > is only 1/7 as big as PROJ, by running du on the source trees with .git > included. > > Pierre > -- > When a barnacle settles down, its brain disintegrates. > Já não percebe nada, já não percebe nada. > > > > ___ > Discuss mailing list > Discuss@lists.osgeo.org > https://lists.osgeo.org/mailman/listinfo/discuss ___ Discuss mailing list Discuss@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/discuss
Re: [OSGeo-Discuss] Parsing well-known text
On Wednesday, 2 September 2020 05:41:40 EDT Even Rouault wrote: > As far as I know, you can't create a CRS WKT with just unit information. The > most minimal content that validates the WKT1 grammar would be something > like: > > LOCAL_CS["unspecified CRS", > LOCAL_DATUM["unspecified datum",2000], > UNIT["metre",1, > AUTHORITY["EPSG","9001"]], > AXIS["Easting",EAST], > AXIS["Northing",NORTH]] That doesn't seem to have an elevation unit; should I assume that they're both meters? > Yes, PROJ >= 6 has support for parsing and creating WKT in several versions > of the WKT standard. > > See proj_create_from_wkt() at > https://proj.org/development/reference/functions.html, and all other > proj_ getters. > > For creation of WKT, you might need the more advanced functions of > https://github.com/OSGeo/PROJ/blob/master/src/proj_experimental.h , before > exporting with proj_as_wkt() That's overkill for what I'm doing. I don't need to make a projection, all I need to know is the units. I don't need to create a WKT except in the case that the original cloud doesn't have one and I'm adding one to specify the units. The code to parse the WKT is in YACC, which I don't know how to work with, and I don't need all the tokens. Treating them as opaque atoms and looking for "UNIT" should work. Where can I find the grammar of the language of which WKT is an instance? Adding any library multiplies the complexity of building and packaging on Windows. PROJ is 79 times as big as Wolkenbase so far and 17 times as big as PerfectTIN. Even Bezitopo, which handles projections, but only conformal ones, is only 1/7 as big as PROJ, by running du on the source trees with .git included. Pierre -- When a barnacle settles down, its brain disintegrates. Já não percebe nada, já não percebe nada. ___ Discuss mailing list Discuss@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/discuss
Re: [OSGeo-Discuss] Parsing well-known text
Pierre, > I just added code to PerfectTIN to read the variable-length records of a LAS > file. I found a WKT record in a square of West Virginia terrain, with no > line feeds; I added line feeds and indentation and attached it. The other > point clouds do not have variable-length records. I'd like to heed the WKT > units when loading a point cloud, if there is a WKT. > > I'm also working on a program called Wolkenbase (not public yet) which will > separate ground from non-ground in a point cloud. If you load a point cloud > with no unit information, I'd like to add a WKT that indicates only the unit > and nothing else. What would that look like? As far as I know, you can't create a CRS WKT with just unit information. The most minimal content that validates the WKT1 grammar would be something like: LOCAL_CS["unspecified CRS", LOCAL_DATUM["unspecified datum",2000], UNIT["metre",1, AUTHORITY["EPSG","9001"]], AXIS["Easting",EAST], AXIS["Northing",NORTH]] > Do you have code to parse and manipulate WKT? The WKT is in a std::string. Yes, PROJ >= 6 has support for parsing and creating WKT in several versions of the WKT standard. See proj_create_from_wkt() at https://proj.org/development/reference/functions.html, and all other proj_ getters. For creation of WKT, you might need the more advanced functions of https://github.com/OSGeo/PROJ/blob/master/src/proj_experimental.h , before exporting with proj_as_wkt() GDAL's OGRSpatialReference class showcases using number of the above mentionned PROJ functions: https://github.com/OSGeo/gdal/blob/master/gdal/ogr/ogrspatialreference.cpp You may also look at PDAL that uses PROJ and/or GDAL underneath to deal with LAS CRS: https://pdal.io/tutorial/las.html#spatial-reference-system Even -- Spatialys - Geospatial professional services http://www.spatialys.com ___ Discuss mailing list Discuss@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/discuss
[OSGeo-Discuss] Parsing well-known text
I just added code to PerfectTIN to read the variable-length records of a LAS file. I found a WKT record in a square of West Virginia terrain, with no line feeds; I added line feeds and indentation and attached it. The other point clouds do not have variable-length records. I'd like to heed the WKT units when loading a point cloud, if there is a WKT. I'm also working on a program called Wolkenbase (not public yet) which will separate ground from non-ground in a point cloud. If you load a point cloud with no unit information, I'd like to add a WKT that indicates only the unit and nothing else. What would that look like? Do you have code to parse and manipulate WKT? The WKT is in a std::string. Pierre -- The Black Garden on the Mountain is not on the Black Mountain. COMPD_CS["NAD83 / UTM zone 17N + NAVD88 height", PROJCS["NAD83 / UTM zone 17N", GEOGCS["NAD83", DATUM["North_American_Datum_1983", SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]], TOWGS84[0,0,0,0,0,0,0], AUTHORITY["EPSG","6269"]], PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]], UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]], AUTHORITY["EPSG","4269"]], PROJECTION["Transverse_Mercator"], PARAMETER["latitude_of_origin",0], PARAMETER["central_meridian",-81], PARAMETER["scale_factor",0.9996], PARAMETER["false_easting",50], PARAMETER["false_northing",0], UNIT["meter",1,AUTHORITY["EPSG","9001"]], AXIS["X",EAST], AXIS["Y",NORTH], AUTHORITY["EPSG","26917"]], VERT_CS["NAVD88 height", VERT_DATUM["North American Vertical Datum 1988",2005,AUTHORITY["EPSG","5103"]], UNIT["meter",1,AUTHORITY["EPSG","9001"]], AXIS["Up",UP], AUTHORITY["EPSG","5703"]]] ___ Discuss mailing list Discuss@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/discuss
