--- a/scripts/compute-homography.py	Thu Apr 02 15:53:59 2015 +0200
+++ b/scripts/compute-homography.py	Wed Apr 08 16:07:15 2015 +0200
@@ -6,8 +6,9 @@
 import numpy as np
 import cv2
 
-import cvutils
-import utils
+import cvutils, utils, storage
+
+# TODO add option to use RANSAC or other robust homography estimation method?
 
 parser = argparse.ArgumentParser(description='The program computes the homography matrix from at least 4 non-colinear point correspondences inputed in the same order in a video frame and a aerial photo/ground map, or from the list of corresponding points in the two planes.', epilog = '''The point correspondence file contains at least 4 non-colinear point coordinates 
 with the following format:
@@ -20,6 +21,7 @@
 in the same order the corresponding points in world and image spaces.''', formatter_class=argparse.RawDescriptionHelpFormatter)
 
 parser.add_argument('-p', dest = 'pointCorrespondencesFilename', help = 'name of the text file containing the point correspondences')
+parser.add_argument('--tsai', dest = 'tsaiCameraFilename', help = 'name of the text file containing the camera parameter following the pinhole camera model (Lund format)') # caution, this is Aliaksei's format
 parser.add_argument('-i', dest = 'videoFrameFilename', help = 'filename of the video frame')
 parser.add_argument('-w', dest = 'worldFilename', help = 'filename of the aerial photo/ground map')
 parser.add_argument('-n', dest = 'nPoints', help = 'number of corresponding points to input', default = 4, type = int)
@@ -78,6 +80,16 @@
 if args.pointCorrespondencesFilename is not None:
     worldPts, videoPts = cvutils.loadPointCorrespondences(args.pointCorrespondencesFilename)
     homography, mask = cv2.findHomography(videoPts, worldPts) # method=0, ransacReprojThreshold=3
+elif args.tsaiCameraFilename is not None: # hack using PDTV
+    f = storage.openCheck(args.tsaiCameraFilename, quitting = True)
+    content = storage.getLines(f)
+    outFilename = '/tmp/camera.yaml'
+    out = storage.openCheck(outFilename, 'w')
+    out.write('data_class: TsaiCamera\n')
+    for l in content:
+        out.write(l.replace(' f:', 'f:').replace(' k:', 'k:').replace(',','.')+'\n')
+    out.close()
+    homography = cvutils.computeHomographyFromPDTV(outFilename)
 elif args.videoFrameFilename is not None and args.worldFilename is not None:
     worldImg = plt.imread(args.worldFilename)
     videoImg = plt.imread(args.videoFrameFilename)
@@ -103,7 +115,7 @@
 if homography.size>0:
     np.savetxt('homography.txt',homography)
 
-if args.displayPoints and args.videoFrameFilename is not None and args.worldFilename is not None and homography.size>0:
+if args.displayPoints and args.videoFrameFilename is not None and args.worldFilename is not None and homography.size>0 and args.tsaiCameraFilename is None:
     worldImg = cv2.imread(args.worldFilename)
     videoImg = cv2.imread(args.videoFrameFilename)
     if args.undistort: