-
Notifications
You must be signed in to change notification settings - Fork 16
/
Copy pathagisoft_all.py
executable file
·155 lines (125 loc) · 4.19 KB
/
agisoft_all.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
#! /usr/bin/env python
#David Shean
#8/22/14
#Script for Agisoft PhotoScanPro workflow
#Based on API v1.0.0, Python 3.3
#Comments include notes for more advanced functionality
#See following forum discussions:
#http://www.agisoft.ru/forum/index.php?topic=2263.0
#http://www.agisoft.ru/forum/index.php?topic=1881.0
import os
import glob
import PhotoScan
import sys
#Need to set the following appropriately
#Path to photos
#photo_fn_path = "/tmp/export"
photo_fn_path = "/Volumes/SHEAN_PHOTO/photo/20140825_MammothTerraces_SfM/export_orig"
photo_fn_ext = "*.jpg"
#Path to ground control file, can contain orientation
gc_fn = "/tmp/gcp.txt"
#Path to calibration file
cal_fn = "/tmp/D800_cal.xml"
#This is the base fn for output files
#out_fn = "/tmp/test"
out_fn = os.path.join(photo_fn_path, "test")
#Define input coordinate system as WGS84
in_crs = PhotoScan.CoordinateSystem()
in_crs.init("EPSG::4326")
#Define output coordinate system as UTM 10N, WGS84
out_crs = PhotoScan.CoordinateSystem()
#out_crs.init("EPSG::32610")
#This is Yellowstone
out_crs.init("EPSG::32612")
#Add timestamp
print("Started")
#Create project file
doc = PhotoScan.app.document
#***
#NOTE: the following (loading photos, ground control, calibration, etc.) can be done manually
#***
#Load photos
new_chunk = PhotoScan.Chunk()
new_chunk.label = "chunk1"
photo_fn_list = glob.glob(os.path.join(photo_fn_path, photo_fn_ext))
for photo_fn in photo_fn_list:
new_chunk.cameras.add(photo_fn)
#Import ground control
gc = new_chunk.ground_control
gc.projection = in_crs
#Load EXIF data from photos
gc.loadExif()
#Alternatively, load csv containing file names and coordinates for photos
#gc.load(gc_fn, "csv")
#Set accuracy of camera positions in meters
#GeoXH
gc.accuracy_cameras = 0.5
#Nikon GP-1
#gc.accuracy_cameras = 5.0
gc.apply()
#Import calibration
#cal = PhotoScan.Calibration(cal_fn)
#new_chunk.calibration_mode('fixed')
#This adds the chunk to the project
doc.chunks.add(new_chunk)
#Update the GUI
PhotoScan.app.update()
doc.save(out_fn + "_init.psz")
#***
#NOTE: end of section with steps that can be accomplished manually
#***
#Grab the active chunk
chunk = doc.activeChunk
#Align photos
print("Aligning photos")
chunk.matchPhotos(accuracy="high", preselection="disabled")
#Use ground control if appropriate for input photos
#chunk.matchPhotos(accuracy="high", preselection="ground control")
chunk.alignPhotos()
PhotoScan.app.update()
doc.save(out_fn + "_sparse.psz")
#NOTE: Adjust bounding box here
#NOTE: markers should be manually identified here
#Build Dense Cloud
print("Building dense cloud")
ncpu = 24
#Not sure about gpu_mask value here, if both cards will be enabled
#Says value 5 enables device number 0 and 2
chunk.buildDenseCloud(quality="medium", filter="mild", gpu_mask=3, cpu_cores_inactive=ncpu)
PhotoScan.app.update()
doc.save(out_fn + "_dense.psz")
#NOTE: Edit dense cloud
#Build Mesh
#NOTE: want to do this both with and without interpolation, export DEM for both
print("Building mesh")
chunk.buildModel(object="arbitrary", source="dense", interpolation="disabled", faces="high")
PhotoScan.app.update()
doc.save(out_fn + "_mesh_nointerp.psz")
#Want to test this smoothing - could help with TIN mesh artifacts
#chunk.smoothModel()
#Build Texture
#chunk.buildTexture(mapping="generic", blending="average", width=2048, height=2048)
#Export DEM
#Should automatically compute appropraite resolution/extent
print("Exporting DEM")
dem_fn = os.path.join(out_fn, "_dem.tif")
chunk.exportDem(dem_fn, format="tif", projeciton=out_crs)
#Export ortho
print("Exporting orthomosaic")
ortho_fn = os.path.join(out_fn, "_ortho.tif")
chunk.exportOrthophoto(ortho_fn, format="tif", blending="average", project=out_crs)
#Export point cloud
print("Exporting point cloud")
#Export las or xyz format
pc_type = "las"
#pc_type = "xyz"
pc_fn = os.path.join(out_fn, "_dense_wgs84.", pc_type)
#Export WGS84 point cloud
chunk.exportPoints(pc_fn, dense=True, precision=7, format=pc_type, projeciton=in_crs)
#Export projected point cloud
#For coord in meters, default precision of 6 is overkill
pc_fn = os.path.join(out_fn, "_dense_proj.", pc_type)
chunk.exportPoints(pc_fn, dense=True, precision=3, format=pc_type, projeciton=out_crs)
#Add timestamp
print("Finshed")