Forum

Author Topic: Compound CRS with geoid works in UI but not in standalone script  (Read 2680 times)

Midkhat

  • Newbie
  • *
  • Posts: 6
    • View Profile
I need to set chunk.crs to 'GDA + AHD', which does exist inbuilt set of CRSs. When using API I need to be able to use geoids from some particular folder, not Metashape one, so I have exact same geoids in Metashape geoid folder and in another one.
When I set CRS in UI (either through Reference Settings menu or with console)  and run Alignment - all works fine.
When I do the same with API - set chunk CRS to 'GDA94 + AHD' and add geoid - MatchPhotos/Alignment fails with the error: "Vertical datum missing". But if I open the file after that - CRS is set to 'GDA94 + AHD' with geoid and Alignment from UI works fine again...
I have checked that the files in Metashape and custom geoids folder are identical.
I have attached two screenshots:
1. Setting CRS in UI - alignment runs fine
2. Opening document after setting CRS from standalone script and failing alignment (if align from this point in UI - works fine)
Code:
Code: [Select]
    crs = [x for x in Metashape.CoordinateSystem().listBuiltinCRS() if x.name == 'GDA94 + AHD height'][0]
    log.debug(f"CRS WKT: {crs.wkt}\n")

    #### get geoid info
    data = gdal.Open(geoid_path)
    metadata = data.GetMetadata()
    wkt = data.GetProjection()
    proj = osr.SpatialReference(wkt=wkt)
    epsg_tiff = proj.GetAttrValue('AUTHORITY', 1)
    log.debug(f"Geoid path: {geoid_path}")
    log.debug(f"Geoid WKT: {wkt}")
    log.debug(f"Geoid EPSG: {epsg_tiff}")
    log.debug(f"Geoid CRS name: {crs.name}")
    log.debug(f"Geoid description: {metadata['TIFFTAG_IMAGEDESCRIPTION']}")
    log.debug(f"CRS WKT: \n{crs.wkt}\n")

    crs.addGeoid(geoid_path)

    chunk.crs = crs
    log.info(f"Chunk <{chunk.label}> CRS is '{chunk.crs.name}'")

Response:
Code: [Select]
2022-07-23 16:59:34,506 - DEBUG - CRS WKT:
COMPD_CS["GDA94 + AHD height",GEOGCS["GDA94",DATUM["Geocentric Datum of Australia 1994",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],TOWGS84[0.06155,-0.01087,-0.04019,0.0394924,0.0327221,0.0328979,-0.009993999999999999],AUTHORITY["EPSG","6283"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9102"]],AUTHORITY["EPSG","4283"]],VERT_CS["AHD height",VERT_DATUM["Australian Height Datum",2005,AUTHORITY["EPSG","5111"]],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","5711"]],AUTHORITY["EPSG","9464"]]

2022-07-23 16:59:34,524 - DEBUG - Geoid path: ./metashape/geoids/au_ga_AUSGeoid09_V1.01.tif
2022-07-23 16:59:34,524 - DEBUG - Geoid WKT: GEOGCRS["GDA94",DATUM["Geocentric Datum of Australia 1994",ELLIPSOID["GRS 1980",6378137,298.257222101,LENGTHUNIT["metre",1]]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],CS[ellipsoidal,3],AXIS["geodetic latitude (Lat)",north,ORDER[1],ANGLEUNIT["degree",0.0174532925199433]],AXIS["geodetic longitude (Lon)",east,ORDER[2],ANGLEUNIT["degree",0.0174532925199433]],AXIS["ellipsoidal height (h)",up,ORDER[3],LENGTHUNIT["metre",1]],USAGE[SCOPE["Geodesy."],AREA["Australia including Lord Howe Island, Macquarie Island, Ashmore and Cartier Islands, Christmas Island, Cocos (Keeling) Islands, Norfolk Island. All onshore and offshore."],BBOX[-60.55,93.41,-8.47,173.34]],ID["EPSG",4939]]
2022-07-23 16:59:34,524 - DEBUG - Geoid EPSG: 4939
2022-07-23 16:59:34,524 - DEBUG - Geoid CRS name: GDA94 + AHD height
2022-07-23 16:59:34,524 - DEBUG - Geoid description: GDA94 (EPSG::4939) to AHD height (EPSG::5711)
2022-07-23 16:59:34,524 - DEBUG - CRS WKT:
COMPD_CS["GDA94 + AHD height",GEOGCS["GDA94",DATUM["Geocentric Datum of Australia 1994",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],TOWGS84[0.06155,-0.01087,-0.04019,0.0394924,0.0327221,0.0328979,-0.009993999999999999],AUTHORITY["EPSG","6283"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9102"]],AUTHORITY["EPSG","4283"]],VERT_CS["AHD height",VERT_DATUM["Australian Height Datum",2005,AUTHORITY["EPSG","5111"]],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","5711"]],AUTHORITY["EPSG","9464"]]

2022-07-23 16:59:34,525 - INFO - Chunk <API_chunk> CRS is 'GDA94 + AHD height'

2022-07-23 16:59:34,933 - CRITICAL - matchPhotos failed.
2022-07-23 16:59:34,933 - CRITICAL - Vertical datum missing



Midkhat

  • Newbie
  • *
  • Posts: 6
    • View Profile
Re: Compound CRS with geoid works in UI but not in standalone script
« Reply #1 on: July 24, 2022, 04:40:21 PM »
I have tried to use the Geoid which is in Applications/Metashape/../geoids folder, same processing fails at Match photos:
Code: [Select]
2022-07-24 23:38:08,579 - DEBUG - Geoid path: /Applications/MetashapePro.app/Contents/MacOS/geoids/au_ga_AUSGeoid09_V1.01.tif
2022-07-24 23:38:08,579 - DEBUG - Geoid EPSG: '4939'
2022-07-24 23:38:08,579 - DEBUG - Geoid description: 'GDA94 (EPSG::4939) to AHD height (EPSG::5711)'
2022-07-24 23:38:08,580 - DEBUG - Chunk CRS set to
2022-07-24 23:38:08,580 - DEBUG - COMPOUND WKT:

COMPD_CS["GDA94 + AHD height",GEOGCS["GDA94",DATUM["Geocentric Datum of Australia 1994",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],TOWGS84[0.06155,-0.01087,-0.04019,0.0394924,0.0327221,0.0328979,-0.009993999999999999],AUTHORITY["EPSG","6283"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9102"]],AUTHORITY["EPSG","4283"]],VERT_CS["AHD height",VERT_DATUM["Australian Height Datum",2005,AUTHORITY["EPSG","5111"]],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","5711"]],AUTHORITY["EPSG","9464"]]

2022-07-24 23:38:08,787 - INFO - Processing step 'MatchPhotos' started.
MatchPhotos: accuracy = Medium, preselection = generic, reference, keypoint limit = 40000, keypoint limit per mpx = 1000, tiepoint limit = 4000, apply masks = 0, filter tie points = 0, filter stationary points = 1, guided matching = 1
2022-07-24 23:38:08,997 - CRITICAL - matchPhotos failed.
2022-07-24 23:38:08,997 - CRITICAL - Vertical datum missing

Alexey Pasumansky

  • Agisoft Technical Support
  • Hero Member
  • *****
  • Posts: 15067
    • View Profile
Re: Compound CRS with geoid works in UI but not in standalone script
« Reply #2 on: July 28, 2022, 12:14:18 PM »
Hello Midkhat,

Does it help, if you add the following line to the beginning of your stand alone script (before any coordinate systems are defined):
Code: [Select]
Metashape.CoordinateSystem.addGeoid("./metashape/geoids/au_ga_AUSGeoid09_V1.01.tif")
addGeoid command should be used to initialize the paths to geoid files that would be later utilized in the script.
Best regards,
Alexey Pasumansky,
Agisoft LLC

Diogo

  • Newbie
  • *
  • Posts: 8
    • View Profile
Re: Compound CRS with geoid works in UI but not in standalone script
« Reply #3 on: December 19, 2022, 06:55:32 PM »
I have the same problem with another crs. Running it in the API from start to beginning, everything works fine. But once I open an aleady processed part the I get the Vertical Daum Missing error. Even when doing the Metashape.CoordinateSystem.addGeoid() command. If I open it via GUI, I can also process the whole project without a problem. Any other sugestions?

Edit
Found out what the problem was:

Metashape.CoordinateSystem.addGeoid() had to be done before opening an older project and after initiating Metashape.document()
« Last Edit: December 20, 2022, 08:15:06 PM by Diogo »