Exporting Chunks as Georeferenced Orthophotos

[AndrewNeubauer]

I am pretty new to python so please forgive me. What I am trying to achieve here is generate orthophotos for all chunks, generating a tif and kmz for chunks which name ends in "RGB" and just a tif for those that end in"NIR". The loop seems to be working correctly as it goes through everyone like it should but there is no file and the console says it exported at 0 x 0 pixels.

My Current Script:

Code: [Select]

from PhotoScan import *

print("This Script will export all chunks into the selected folder")

doc = PhotoScan.app.document


CoordinateSystemEPSG = "EPSG::4326"
CoordinateSystem = PhotoScan.CoordinateSystem()
CoordinateSystem.init(CoordinateSystemEPSG)

def getDir(name):
dir = PhotoScan.app.getExistingDirectory(name)
if dir == '':
print("Must provide a valid directory")
return
return dir


saveFolder = getDir("Select Save Folder")




print("---Exporting Orthoimage(s)...")

for chunk in doc.chunks:
file = saveFolder + "\\"  + chunk.label
print("File: " + file)

if chunk.label[-3:] == "RGB":
chunk.exportOrthophoto(file + ".tif", format="tif", blending="mosaic", projection=CoordinateSystem, dx=0.2, dy=0.2)
chunk.exportOrthophoto(file + ".kmz", format="kmz", blending="mosaic", projection=CoordinateSystem, dx=0.05, dy=0.05, write_kml= True)
elif chunk.label[-3:] == "NIR":
chunk.exportOrthophoto(file + ".tif", format="tif", blending="mosaic", projection=CoordinateSystem, dx=0.2, dy=0.2)
else:
app.messageBox("Exporting orthoimage failed: " + chunk.label)
The console spits out this (multiple times):

Code: [Select]

File: ****\** RGB
initializing renderer... tessellating mesh...done (335152 -> 335152 faces)
done in 0.338 sec
Raster size: 0x0
Finished processing in 0.355 sec (exit code 0)
initializing renderer... tessellating mesh...done (335152 -> 335152 faces)
done in 0.364 sec
Raster size: 0x0
Finished processing in 0.364 sec (exit code 0)
File: ****\** NIR
initializing renderer... tessellating mesh...done (342265 -> 342265 faces)
done in 0.295 sec
Raster size: 0x0
Finished processing in 0.295 sec (exit code 0)

[Alexey Pasumansky]

Hello AndrewNeubauer,

For the exportOrthophoto function you need to input the resolution in the units used for the current projection. So for WGS84 you need to specify the resolution in degrees, not in meters. You can take the suggested values from Export Orthophoto dialog of GUI.

[AndrewNeubauer]

Thanks for the reply Alexey,

I did not realize that, I need to export using resolution in meters, and all of the chunk sizes are different. How should I go about doing that? Am I going to have to calculate that in the script, and if so, how?

[Alexey Pasumansky]

Hello Andrew,

If chunks refer to almost the same area and all are referenced in WGS84, you can take x and y resolution from Export Orthophoto dialog (in degrees) for any sample project for the same area.

To calculate the resolution using Python I can suggest the following method:

Code: [Select]

#recalculating WGS84 resolution from degrees into meters
crd = chunk.ground_control.locations[chunk.cameras[0]].coord  #first image coordinates, should be aligned

#longitude
v1 = PhotoScan.Vector((crd[0], crd[1], 0) )
v2 = PhotoScan.Vector((crd[0] + 0.001, crd[1], 0))
vm1 = chunk.crs.unproject(v1)
vm2 = chunk.crs.unproject(v2)
res_x = (vm2 - vm1).norm() * 1000

#latitude
v2 = PhotoScan.Vector( (crd[0], crd[1] + 0.001, 0))
vm2 = chunk.crs.unproject(v2)
res_y = (vm2 - vm1).norm() * 1000

pixel_x = pixel_y = 0.2  #export resolution (meters/pix)
d_x = pixel_x / res_x 
d_y = pixel_y / res_y

[AndrewNeubauer]

Alexey,

This may not be the most pythonic or efficient way to do this but it is working great for me, I put it into a function to make this super simple. Thank you so much!

Code: [Select]

def toMeters(res):
if res != float(res):
return 0
else:
try:
#recalculating WGS84 resolution from degrees into meters
crd = chunk.ground_control.locations[chunk.cameras[0]].coord  #first image coordinates, should be aligned

#longitude
v1 = PhotoScan.Vector((crd[0], crd[1], 0) )
v2 = PhotoScan.Vector((crd[0] + 0.001, crd[1], 0))
vm1 = chunk.crs.unproject(v1)
vm2 = chunk.crs.unproject(v2)
res_x = (vm2 - vm1).norm() * 1000

#latitude
v2 = PhotoScan.Vector( (crd[0], crd[1] + 0.001, 0))
vm2 = chunk.crs.unproject(v2)
res_y = (vm2 - vm1).norm() * 1000

#export resolution (meters/pix)
d_x = res / res_x 
d_y = res / res_y
return [d_x, d_y]

except IndexError:
return 0Returning 0 means it failed, otherwise it returns the x and y resolution in a list.

[marco]

Hi,
I tried this solution too, but I always get the following error:

crd = chunk.ground_control.locations[chunk.cameras[0]].coord
TypeError: invalid key

If I output the variables everything seems to be ok.

Code: [Select]

    print(chunk.cameras[0])
    print(chunk.ground_control.locations)

    crd = chunk.ground_control.locations[chunk.cameras[0]].coord
    [...]
<Camera 'test0001.JPG'>
<PhotoScan.GroundControlLocations object at 0x7ff9a416e1f0>

Any obviously mistake? Or is there a special precondition for accessing the ground_control.locations?

Thanks in advance

[Alexey Pasumansky]

Hello marco,

chunk.ground_control.locations contain only those items, that have source values input in the Ground Control pane.

[giancan]

Dear all,
I tried the suggested code or exporting orthophoto on my project (in metric coordinate system) and the processing seems fine, included the info
Raster size: 12975x14829

The only problem is that I don't get any file!!!!!
I tried to export in several locations (in case it was a path problem) but nothing changed.
Any suggestion?
Thanks

[Alexey Pasumansky]

Hello giancan,

Could you please post the export orthophoto line from your code? If it uses the variable for the path, please also include the variable value. Probably, the path or filename is incorrect.