opengeospatial / geotiff

18 stars 9 forks source link

Handling projected CRS + ellipsoidal height #76

Closed EmDevys closed 4 years ago

EmDevys commented 5 years ago

From Max Martinez (email 30/03) to GeoTIFF.SWG D.3. Geographic 3D CRS (case of ellipsoïdal height) (page 84) The recommendation as currently worded does not seem to encompass projected + ellipsoidal height. Proposed re-wording: a) give GTModelTypeGeoKey value of 1 or 2 and imply or identify a geographic 2D CRS, giving the geographic 3D CRS description through the VerticalGeoKey. The implied or specified geographic 2D CRS should be the horizontal component of the geographic 3D CRS. (This is the recommended option)

From Roger Lott (01/04) We have discussed this. My recollection is that what Max describes would be considered in v2.0 to consider, as part of a more general consideration of CRS subtypes not covered in v1.0. (In one draft I suggested some others that are missing, e.g. compound of various forms). So for v1.1 I think the text is ok as on P84.

From Max (01/04) I don’t recall there being any presentation of evidence of use of compound in existing GeoTIFF community practice. Thus the decision to leave that to v2. Whereas DGIWG and other standards explicitly endorse the use of 4979 in contexts where either PCS or GCS can be specified (at least I couldn’t find any mention of use of 4979 precluding the use of a PCS in that spec….maybe I missed it). Once v1.1 assents to describing ellipsoidal height via a 3D geographic CRS in the vertical key (even if only in annex), recommending that only when the horizontal is explicitly 2d geographic (and not in the case of the derived PCS) seems pointless. We’re not protecting anybody from surprises by doing so yet we’re making people who wish to do so unsure of its validity.

From Roger (01/04) OGC Topic 2 = ISO 19111 explicitly states that an ellipsoidal height cannot exist alone but only as a part of a 3D system. So ellipsoidal height cannot be considered as part of a compound CRS. It was for this reason that the projected 3D CRS was added to the recent Topic 2 revision. When I suggested that this CRS subtype be added to v1.1 (with some other subtypes also missing) they were considered out of scope of this minor revision and to be considered for v2.0. Those other subtypes included geog3D CRS. The current provisions are a work-around these decisions and allow (i) retention of the 5 values allowed for GTModelTypeKeoKey for CRS type (unknown, geog2D, geocentric, projected [2D] or user-defined as given in v1.0, one of which is mandatory, (ii) no extension of this list until v2.0, and (iii) In v1.1 the use Max describes through citing the geog3D CRS code in the Vertical GeoKey. The name 'VerticalGeoKey' was deliberately selected so that it omitted any CRS implication [all other analogous keys include CRS in their name).. The appropriate place to be identifying the type of CRS is through the GTModelTypeKeoKey. We have wilfully decided that extending the v1.0 list was out of scope of this revision. Not what I would have decided if I had unilaterally written v1.1, but I am happy to go along with the majority decision that changes of this nature are out of scope of the v1.1 revision. So I think we leave this part of the v1.1 document as it is.

If at this rather late stage we do decide to reverse our earlier decision on this, then I feel that the extension should not just be for Geog3D but for all other missing types that I listed some time ago. This is not a one or two sentence change, but will require some substantive supporting text in a couple of areas of the document.

Max (02/04) Seems like a strange result. Note that I never suggested extending the allowed values for the GTModelTypeGeoKey. There seemed no need to reverse any previous decision in order to tweak the wording in D.3 so that it applies whenever there is a 2D geographic specified or implied in the Model CRS.

A picture might help here. Below is the representation of how the EPSG database seems to relate CRS definitions to their BaseCRS.

ModelType == 1: 2D Projected  2D Geographic  3D Geographic  3D Geocentric ModelType == 2: 2D Geographic - 3D Geographic  3D Geocentric

D3 says you can properly indicate ellipsoidal height by putting the highlighted CRS code in the VerticalGeoKey…but only if the ModelType is 2.

EmDevys commented 5 years ago

Could we decide on whether we leave this issue opened till v2 (and keep as is in v1.1) ? Or bring some additional note (or a BP) allowing the use of VerticalGeoKey with a Geographic 3D CRS, for 22 projected as well, which seems a possible use case (e.g allowed by DGIWG on basis of GeoTIFF 1.0). If not, 1.1 would be more strict than 1.0. In addition, there is no requirement in 1.1 excluding this, just the informative recommendations provided in Annex D. To be further discussed

EmDevys commented 4 years ago

Roger Lott (email 25 June) DGIWG's comment asks for a ellipsoidal height in a projected CRS. In the revised 2019 Topic 2 model, this would be a projected 3D CRS. That construct was not available to GeoTIFF v1.0 and it is therefore out of scope for v1.1.

DGIWG is correct in saying that "nothing excludes the use of an ellipsoidal Vertical CS Type Code for VerticalCSTypeGeoKey with a GeographicTypeGeoKey or a ProjectedCSTypeGeoKey". However GeoTIFF v1.0 in §6.3.4.1 lists numerous non-existent EPSG codes. So nothing in the format excludes the use of a code, but no codes exist.

The underlying flaw in the v1.0 document is assuming that a 1D vertical CRS may be either gravity-related or ellipsoidal. In the Topic 2 data model ellipsoidal height can only be an axis in a 3D geographic CRS, as explained in a couple of places in v1.1 and in particular D.3. A height can be associated with a projected CRS as long as it is a gravity-related height (then it is described implicitly as a compound CRS, projected CRS with a vertical CRS). This is clearly explained in the v1.1 document.

In conclusion, support for ellipsoidal height with projected CRS is out of scope for v1.1 so has to await v1.2 or later.

EmDevys commented 4 years ago

Roger Lott (email 26 June) Unlike some of the geographic and projected CRS codes given in v1.0 which were at some time published by EPSG and then withdrawn because errors were found in them, so now invalid (deprecated), the vertical ellipsoid height codes listed in GeoTIFF v1.0 never existed in EPSG. There never were EPSG codes – let alone valid codes. You cannot invalidate something that does not exist. I did meet with Ritter and Ruth in DC to discuss the inclusion of EPSG codes in GeoTIFF and this was never discussed. I can't explain how they got there. Had the GeoTIFF v1.0 authors approached EPSG to create records for these codes, we would never have agreed to it. The concept is not supported in the 19111/Topic 2 data model, and never has been, nor in the EPSG data model that preceded the introduction of the ISO model. It is not something recognised by the geodetic community.

PS Not directly related to GeoTIFF, and certainly not to v1.1, I would like to take the opportunity to discuss with the RS community what exactly is done in imagery processing to create products that would use a projected 3D (easting, northing, ellipsoidal height) CRS. How is the ellipsoidal height value generated away from the sensor and image principle point/nadir? The issue impacts whether, and if so how, EPSG might describe projected 3D CRSs.

EmDevys commented 4 years ago

Backward compatibility issue added to Annex H. The SWG has taken the action to provide an update: Backward compatibility annex should refer to Table D.1. Revision 1.1 does not support ellipsoidal height with projected CRS (because the GeoTIFF mechanisms in GeoTIFF 1.0 and 1.1 do not allow any combined CRS mechanism and there is no code for ellipsoidal height in EPSG register). => potential future work item for 1.2

max-martinez commented 4 years ago

Roger (apparently) wrote:

In the revised 2019 Topic 2 model, this would be a projected 3D CRS. That construct was not available to GeoTIFF v1.0 and it is therefore out of scope for v1.1.

You only need a projected 3D CRS to support a 3D grid. I think it needs to be explained under what conditions you would be interpreting a raster DEM in TIFF file as a 3D Grid.

I would interpret it as a 2D Grid with an attribute that represents height, the interpretation of which is guided by metadata in the VerticalCRSGeoKey.

Annex D.2 notwithstanding, I suspect that close to 100% of the GeoTIFF accessing implementations out there would hope to never find a Compound CRS indicated in the GTModelTypeGeoKey of the GeoTIFF (or a 3D Geographic CRS present in the GeographicTypeGeoKey, for that matter). All that represents is a nuisance. These implementations want to access a 2D Grid and, if interested in raster DEMs, have the necessary metadata available to properly interpret the heights represented by the pixel values. Whether the GTModelTypeeGeoKey is Projected or Geographic, this can be achieved by putting either a VerticalCRS or 3D GeogrpahicCRS EPSG code in the VerticalCRSGeoKey.

Seeing as this is:

  1. Not prohibited by the GeoTIFF 1.1 specification (though not clearly advised in D3 due to our inaction)
  2. Seems to be consistent with the EPSG database itself as suggested in an earlier comment in this thread
  3. Seems to be an already existing practice of DGIWG
  4. Has been resisted only on the grounds that it requires a projected 3D CRS (which it doesn't)

Hexagon will be moving its software implementation toward both interpreting and populating (when so directed) GeoTIFFs with 2D projected GTModelTypeGeoKeys and 3D Geographic VerticalCRSGeoKeys as projected raster DEMs whose heights are expressed as the relevant ellipsoidal height.

EmDevys commented 4 years ago

Dear Max Thanks for comment.

This has been a strong topic of discussion, and it was decided to keep Annex D.3 (Informative) as is. As you remind it, there is nothing prohibiting the use of a VerticalGeoKey with any kind of GTModelTypeGeoKey in the GeoTIFF 1.1 requirements.

But Annex D.3 provides restrictions as indicated in the present GeoTIFF 1.1.

As you remind it, DGIWG GeoTIFF profile allows the EPSG code 4979 (Geographic 3D WGS 84 ellipsoid) whatever the GTModelTypeGeoKey. However, EPSG registry specifies that the Coordinate System with this code is Geographic 3D, lat, long, ellipsoidal height. Therefore this does not seems to be an acceptable solution, and DGIWG has to confirm such a use case.

However, there might be use cases with ellipsoidal height and projected CRS, so there should be some best practice for addressing this, by use of VerticalGeoKey and/or VerticalDatumGeoKey.

I copy the mailing-list for information and awareness. Any comment is welcome.

By the way, as the ballot deadline is today, you may wish to send this comment as a formal ballot comment. Please check with Frederic who sent the Hexagon vote (I copy Frederic to this email as well), which is presently a yes vote without comment.

Thanks in advance

Emmanuel

De : max-martinez [mailto:notifications@github.com] Envoyé : vendredi 23 août 2019 23:10 À : opengeospatial/geotiff Cc : Emmanuel Devys; Author Objet : Re: [opengeospatial/geotiff] Handling projected CRS + ellipsoidal height (#76)

Roger (apparently) wrote:

In the revised 2019 Topic 2 model, this would be a projected 3D CRS. That construct was not available to GeoTIFF v1.0 and it is therefore out of scope for v1.1.

You only need a projected 3D CRS to support a 3D grid. I think it needs to be explained under what conditions you would be interpreting a raster DEM in TIFF file as a 3D Grid.

I would interpret it as a 2D Grid with an attribute that represents height, the interpretation of which is guided by metadata in the VerticalCRSGeoKey.

Annex D.2 notwithstanding, I suspect that close to 100% of the GeoTIFF accessing implementations out there would hope to never find a Compound CRS (or a 3D Geographic CRS, for that matter) in the GTModelTypeGeoKey of the GeoTIFF. All that represents is a nuisance. These implementations want to access a 2D Grid and, if interested in raster DEMs, have the necessary metadata available to properly interpret the heights represented by the pixel values. Whether the GTModelTypeeGeoKey is Projected or Geographic, this can be achieved by putting either a VerticalCRS or 3D GeogrpahicCRS EPSG code in the VerticalCRSGeoKey.

Seeing as this is:

  1. Not prohibited by the GeoTIFF 1.1 specification (though not clearly advised in D3 due to our inaction)
  2. Seems to be consistent with the EPSG database itself as suggested in an earlier comment in this thread
  3. Seems to be an already existing practice of DGIWG
  4. Has been resisted only on the grounds that it requires a projected 3D CRS (which it doesn't)

Hexagon will be moving its software implementation toward both interpreting and populating (when so directed) GeoTIFFs with 2D projected GTModelTypeGeoKeys and 3D Geographic VerticalCRSGeoKeys as projected raster DEMs whose heights are expressed as the relevant ellipsoidal height.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHubhttps://github.com/opengeospatial/geotiff/issues/76?email_source=notifications&email_token=ACDHG2VIW3FNM56GIHCZ3CTQGBGZ7A5CNFSM4HC563P2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD5BKZPA#issuecomment-524463292, or mute the threadhttps://github.com/notifications/unsubscribe-auth/ACDHG2UQRSSNBOW2EEXJRTLQGBGZ7ANCNFSM4HC563PQ.

max-martinez commented 4 years ago

As you remind it, DGIWG GeoTIFF profile allows the EPSG code 4979 (Geographic 3D WGS 84 ellipsoid) whatever the GTModelTypeGeoKey. However, EPSG registry specifies that the Coordinate System with this code is Geographic 3D, lat, long, ellipsoidal height. Therefore this does not seems to be an acceptable solution, and DGIWG has to confirm such a use case.

This seems problematic only if one insists that the information in the VerticalGeoKey must perform some function in the definition of the model space of the raster. Although this sounds like a reasonable stance when the model space is 3D, it is a non-issue when the model space is 2D. This is and has been one of the sources of disagreement here. Many SWG members seem to be champing at the bit to get 3D CRS definitions into the model type (Compound, 3D Geographic, and, if/when AS 2/ISO 19111 allows, projected with ellipsoidal). I thank them for their restraint in maintaining backward compatibilty with 1.0. Moving beyond 1.1, however, I view the introduction of 3D CRS definitions in other than the VerticalGeoKey (via extensions to the GTModelTypeGeoKey) for support of raster DEMs as unnecessary and ultimately detrimental to interoperability.

In our use cases for raster DEMs, we view the VerticalGeoKey as metadata for the attribute value of the 2D coverage that is the GeoTIFF dataset. As noted in the past, SWE Common (see referenceFrame, axisID in section 7.2.3 describing AbstractSimpleComponent) appears to allow this sort of per-attriubte metadata not only for raster dem values but for any coordinate value of any CRS one can devise.

Most of our raster DEM use cases boil down to, "given a ground X, Y, what is the elevation at that location?" We don't start with an I, J. And we don't have a Z. How is a 3D CoordinateOperation going to help us here? With a 2D CoordinateOperation and the VerticalGeoKey metadata, we can ultimately get to where we need to go. Whether the coordinate operations we need to apply are conversions or, possibly to the dismay of geodisists, transformations, we would like to leave those decisions to our users.

I did notice that in Annex F of GeoTIFF 1.1 the original GeoTIFF 1.0 example with Z scaling is illustrated and an additional DGIWG example also is making use of the 3D model space parameters. Thus, one way to be as flexible as possible is to give GeoTIFFs with a Sz = 1.0 and Tz = 0.0 (or the transformation matrix equivalent) the same treatment as 2D Grids with a VerticalGeoKey. Likely the only way to make Z scaling that is not a no-op work properly for a raster access API that returns a 2D->2D georeferencing transformation is to maintain the definition found in the VerticalGeoKey (to be returned as metadata) but to apply the scale before returning pixel values from the data access API. This will likely require synthesizing a pixel type that is not a direct translation of the BitsPerSample+SampleFormat information in some cases so as to be able to represent the scaled data without loss. It would have been nice to not have to support this at all since people seem to have gotten along without it to this point, but seems Even has added this to GDAL so likely software developers will have to put something in to handle it now that the community will be able to generate these types of GeoTIFFs.

EmDevys commented 4 years ago

Max thanks for additional clarification on this view of the GeoTIFF VerticalCSTypeGeoKey (in 1.0), renamed VerticalGeoKey in proposed 1.1, as "metadata for the attribute value of the 2D coverage that is the GeoTIFF dataset" in case of elevation/depth data. The problem I identify (if I may) is that we (in GeoTIFF) are using EPSG CRS codes (such as EPSG code 3855 or 4979 for this), and therefore objects / terms / tools provided by geodesists. We (this SWG) should be careful to use terms and state the requirements correctly (without ambiguity no potential conflict in interpretation). But I must confess that I personnally feel close to your opinion and proposal in principle.

Now the key question is: how to solve this disagreement and finalize GeoTIFF 1.1?

(**): I also note your point that "the introduction of 3D CRS definitions in other than the VerticalGeoKey (via extensions to the GTModelTypeGeoKey) for support of raster DEMs as unnecessary and ultimately detrimental to interoperability." This view is a keypoint and key to the evolution of GeoTIFF. @SWG: please feel free to provide your view. We must finalize the decision on this before publication (so hopefully by Banff OGC TC).

Z scaling-factor (3rd component of ModelPixelScaleTag): seems to be another topic of discussion

EmDevys commented 4 years ago

Roger Lott (3/09/2019) To meet Max's suggestion I think D.3 needs to be changed to include a new b) inserted after a) :-

• b) give GTModelTypeGeoKey value of 1 and identify a projected CRS through the ProjectedCRSGeoKey, giving the geographic 3D CRS description through the VerticalGeoKey. The projected CRS should have as its base CRS the geographic 2D CRS which is the horizontal component of the geographic 3D CRS.

• Then change current b) to c)

• And in the paragraph before the heading to table D.1, change 'Neither of these …' to 'None of these …'

From a geodetic viewpoint this is no different (or, depending upon your perspective, just as bad) as allowing a). A measurement is not a coordinate.

I don't think 7.4.4 needs any change.

For both a) and this proposed new b) there can be cases where the geographic 2D is not derived, so there is no geographic 3D – ED50, NAD27 and NAD83 are examples. In the current EPSG Dataset there are 544 geographic 2D CRSs but only 181 of these are derived from a geographic 3D. The only saving grace is that those that are derived from a geographic 3D are the modern CRSs including WGS 84, so more likely to be used than the older CRSs.

EmDevys commented 4 years ago

Emmanuel Devys (3/09/2019) Thanks for this proposal, which seems to me an acceptable compromise on this issue. Hoping everybody agrees on this, if so, I’ll send out to the SWG (and also add it under the github for traceability) for information + validation (under the usual: is there any objection …), before adjusting the draft.

For the Annex D, you submit the question requirement vs. recommendation. Presently in the document, all requirements are in Chapter 7, and the title of Annex D is : recommendations (there is even an “informative” into brackets). In D.2, the technical statement is preceded by “Recommendation” (which is not the case in D.3). I agree on your worry on “what happens if the condition stated” is not respected? (May be there should be an additional warning NOT to do this ? (or any more solution that you / we believe would be more appropriate).

EmDevys commented 4 years ago

Max Martinez (3/09/2019) Summary: I'm OK with Roger's proposed revision to D.3. I don't think the warning is necessary, but I consent if it makes the proposal more palatable to others.

Cath-up details: I had actually proposed modified wording to D.3 back in March (see https://github.com/opengeospatial/geotiff/issues/76?_pjax=%23js-repo-pjax-container#issue-428101889). As far as I am concerned, this is equivalent to what Roger has just proposed so I am fine with his wording.

I feel 'should' is the right term. This is a best practice we are describing in an annex, not a requirement we are asserting in the specification.

So, yes, if someone happens to use or imply a 2D geographic for the horizontal coordinates that is unrelated to the 3D geographic CRS specified in the VerticalGeoKey, then a transformation would have to be performed to produce the Lat, Lon that represent the horizontal coordinate values that are to be combined with the Z to get the 3D geographic coordinate. That would introduce some error if that were required by one’s data arrangement. So you shouldn't do it.

Roger wrote: We will just leave it to software to throw exceptions if the two CRSs do not match

In software development, we try to recognize the difference between "unworkable" and "inadvisable" so as to be as flexible as possible in our implementations. This is not an "unworkable" situation for us, only "inadvisable" so we would be inclined to issue a warning, not throw an exception.

Roger wrote: It seems to me that we may have a problem with common understanding of whether a height value is a coordinate or a measurement. CRSs give definition to coordinates. If height is referenced in a CRS definition, it is considered to be a coordinate.

Yes, we might not have a common understanding here. I think it would be worth making the observation at this point that a dataset containing merely 3D coordinate values is of zero value to anyone. Do you agree?

Roger wrote: But is it actually a measurement being used as an attribute?

In the absence of specified semantics that assert otherwise, I interpret a raster DEM to be an assertion that the ground elevation at each 2D grid coordinate has been measured to be the range value that corresponds to that grid point. The range value is to be interpreted in the context of the provided metadata (in this case, the VerticalGeoKey value). So, yes, I believe it is measurement being used as an attribute.

Roger wrote: Now, are you going to tell me that latitude and longitude (or easting and northing) are also considered to be attributes of a cell ??

Er, no. Attributes are found in the range. The lat,lon (or E,N) are related to the domain of the 2D grid coverage. For a rectified grid, this relation is via the coordConv (please see AS Topic 6 section 8.9.4).

Roger wrote: AS Topic 2 is quite clear that ellipsoidal height cannot exist in isolation, only as part of a 3D CRS. If that is thought to be wrong then we need to justify why, and have Topic 2 changed.

I have no issue with this. What I am unclear about is why you think that applies to this case (assuming we are talking about a georectified dataset). A georectified raster has grid points whose horizontal geospatial location are fully known. The Z measurements at those grid points can be combined with the appropriate Lat, Lon coordinates computed from those locations to establish an acceptable 3D geographic coordinate, no?

Roger wrote: I am intrigued to understand how ellipsoidal height is obtained when it is associated with projected 2D coordinates.

I believe I just explained it above, but, again, it is worth pointing out that in every case we can consider here (projected + ellipsoidal, geographic + ellipsoidal, projected + vertical, geographic + vertical), there are no explicit horizontal ground coordinates. There aren't even explicit grid coordinates. That's the secret sauce of a raster...the domain is implicit. So no matter which case we consider, we need to perform some sort of coordinate conversion to produce horizontal coordinates that can be married to the Z values to give a 3D ground coordinate. What is unclear to me is why a coordinate conversion that involves a projection is more offensive than one that does not. I thought AS Topic 2/ISO 19111 ultimately abstracted all lossless coordinate operations under the same umbrella.

With respect to z-scaling, I have sympathy with the original specification authors who, in 1995 (and earlier) were concerned with being able to be as efficient as possible in terms of the size of datasets. That's why I assume they put in the provision to allow well communicated information on scaled values (that would allow the use of a more space efficient bits-per-sample under some circumstances). I was merely observing that I don't believe this aspect of the specification to have been widely used and it might require some awkward decisions for software developers to serve this data under their current APIs. I think the provision is clear as originally specified and we just need to make those awkward decisions for our implementations (either fudge the CRS, fudge the pixel type and pixel values, or leave it to the users of the API to resolve the disconnect by using additional metadata that they may not have been using before).

cmheazel commented 4 years ago

Closed through consent of the GeoTIFF membership at the Banff TC meeting