end harris

TPs/harris.cpp

@@ -5,6 +5,7 @@
 #include <fstream>
 #include <vector>
 #include <algorithm>
+#include <assert.h>
 
 using namespace std;
 using namespace Imagine;
@@ -12,6 +13,7 @@ using namespace Imagine;
 /************************/
 /****** PARAMETERS ******/
 /************************/
+
 const static double focal = 2600;
 
 // Standart deviation of the gaussian used for the gradient
@@ -32,7 +34,10 @@ const static double rNMS = 5;
 const static double c = 0.9;
 
 // Number of points to plot
-const static int N = 200;
+const static int N = 800;
+
+// Number of points to be considered for location error
+const static int NLEs[3] = {100, 300, 600};
 /************************/
 
 // Return a 1D gaussian kernel of standard deviation sigma of size 2k+1
@@ -242,8 +247,8 @@ FMatrix<double,3,3> getScale(double scale) {
 // Return a random homography defined by a rotation of the camera
 FMatrix<double,3,3> getRandomHomography(double w, double h) {
 	FMatrix<double,3,3> Rx, Ry, Rz, K, H;
-	Rx = getRotation(1.0 * (doubleRandom() - 0.5));
-	Ry = getRotation(1.0 * (doubleRandom() - 0.5));
+	Rx = getRotation(1.1 * (doubleRandom() - 0.5));
+	Ry = getRotation(1.1 * (doubleRandom() - 0.5));
 	Rz = getRotation(0.3 * (doubleRandom() - 0.5));
 	swap(Ry(0, 1), Ry(0, 2));
 	swap(Ry(1, 0), Ry(2, 0));
@@ -346,12 +351,23 @@ vector<pair<double, int>> ANMS(const vector<pair<double, int>> &points, int w) {
 	return res;
 }
 
-// Plot the N best points
-void plotPoints(const vector<pair<double, int>> &points, int N, int w, double fact=1.0, int ox=0, int oy=0) {
+// Plot the N best points and also save image
+void plotPoints(const vector<pair<double, int>> &points, int N, Image<byte> Ib, string name, double fact=1.0, int ox=0, int oy=0) {
+	int w = Ib.width(), h = Ib.height();
+	Image<Color> Ip(w, h);
+	#pragma omp parallel for
+	for(int i = 0; i < w; i++) 
+		for(int j = 0; j < h; j++)
+			Ip(i, j) = Color(Ib(i, j));
 	for(int i = 0; i < N; i++) {
 		int p = points[i].second;
-		drawCircle(ox + fact*(p % w), oy + fact*(p / w), 4, RED, 3);
+		int x = p % w, y = p / w;
+		drawCircle(ox + fact*x, oy + fact*y, 4, RED, 3);
+		for(int a = max(0, x-3); a <= min(w-1, x+3); a++)
+			for(int b = max(0, y-3); b <= min(h-1, y+3); b++)
+				Ip(a, b) = Color(255, 0, 0);
 	}
+	save(Ip, name);
 }
 
 // Save the responses in res in a file 'name' to plot the number of detection in funcion of the threshold
@@ -367,11 +383,15 @@ bool save_responses(const vector<pair<double, int>> &res, string title) {
 	return true;
 }
 
+// Compute the minimal distance for the N best points of res2 multiplied by H^{-1} to the N bets points of res
 bool location_error(const vector<pair<double, int>> &res, vector<pair<double, int>> &res2,
-					double x0, double y0, FMatrix<double,3,3> H, int w, int w2, string title) {
+					double x0, double y0, FMatrix<double,3,3> H, int w, int w2, string title, int N)
+{
+	assert(N <= int(res.size()) && N <= int(res2.size()));
+	
 	// We sort points in function of the abcisse
 	vector<pair<double, double>> points;
-	for(auto &a : res) points.push_back({a.second % w, a.second / w});
+	for(int i = 0; i < N; i++) points.push_back({res[i].second % w, res[i].second / w});
 	sort(points.begin(), points.end());
 	int n = points.size();
 
@@ -379,8 +399,8 @@ bool location_error(const vector<pair<double, int>> &res, vector<pair<double, in
 	v[2] = 1;
 	vector<double> ts;
 	FMatrix<double,3,3> invH = inverse(H);
-	for(auto p : res2) {
-		v[0] = p.second % w2 + x0, v[1] = p.second / w2 + y0;
+	for(int i = 0; i < N; i++) {
+		v[0] = res2[i].second % w2 + x0, v[1] = res2[i].second / w2 + y0;
 		v = invH * v;
 		v /= v[2];
 		double x = v[0], y = v[1];
@@ -417,16 +437,15 @@ bool location_error(const vector<pair<double, int>> &res, vector<pair<double, in
 				}
 			}
 		} while(!stop);
-		ts.push_back(dist);
+		ts.push_back(sqrt(dist));
 	}
 	sort(ts.begin(), ts.end());
 
 	// Then save
-	int n2 = res2.size();
 	ofstream f;
-	f.open(title + "_loc.txt");
+	f.open(title + "_" + to_string(N) + "_loc.txt");
 	if(!f.is_open()) return false;
-	f << title << "\n" << n2 << "\n";
+	f << title << "\n" << N << "\n";
 	for(double t : ts) f << t << "\n";
 	f.close();
 	return true;
@@ -441,33 +460,40 @@ void compute_NMS_and_display(Window W, int i, Image<byte> Ib, vector<pair<double
 {
 	setActiveWindow(W, i);
 	display(Ib, ox, oy, false, fact);
-	plotPoints(responses, N, w, fact, ox, oy);
+	plotPoints(responses, N, Ib, name+".jpeg", fact, ox, oy);
 	save_responses(responses, name);
 	applyThreshold(responses);
 	if(res0.empty()) res0 = responses;
-	else location_error(res0, responses, x0, y0, H, w0, w, name);
+	else for(int NLE : NLEs)
+		location_error(res0, responses, x0, y0, H, w0, w, name, NLE);
 
 	// best NMS responses
 	setActiveWindow(W, i+5);
 	display(Ib, ox, oy, false, fact);
 	vector<pair<double, int>> nms = NMS(responses, w, h, R, offset);
-	plotPoints(nms, N, w, fact, ox, oy);
+	plotPoints(nms, N, Ib, name+"_nms.jpeg", fact, ox, oy);
 	if(nms0.empty()) nms0 = nms;
-	else location_error(nms0, nms, x0, y0, H, w0, w, name+"_nms");
+	else  for(int NLE : NLEs)
+		location_error(nms0, nms, x0, y0, H, w0, w, name+"_nms", NLE);
 
 	// best ANMS responses
 	setActiveWindow(W, i+10);
 	display(Ib, ox, oy, false, fact);
 	vector<pair<double, int>> anms = ANMS(responses, w);
-	plotPoints(anms, N, w, fact, ox, oy);
+	plotPoints(anms, N, Ib, name+"_anms.jpeg", fact, ox, oy);
 	if(anms0.empty()) anms0 = anms;
-	else location_error(anms0, anms, x0, y0, H, w0, w, name+"_anms");
+	else  for(int NLE : NLEs)
+		location_error(anms0, anms, x0, y0, H, w0, w, name+"_anms", NLE);
 }
 
-int main() {
+int main(int argc, char* argv[]) {
+	if(argc != 2) {
+		cerr << "Exactly one argument needed, corresponding to the path of the image" << endl;
+		return 1;
+	}
 	initRandom();
 	Image<byte> I0;
-	if(! load(I0, srcPath("im1.jpg"))) {
+	if(! load(I0, argv[1])) {
 		cerr << "Unable to load the image !" << endl;
 		return 1;
 	}
@@ -521,12 +547,12 @@ int main() {
 	// best Harris responses on noised image
 	cout << "Noise" << endl;
 	H = FMatrix<double,3,3>::Identity();
-	double sigma = 20;
+	double sigma = 60;
 	I2 = applyNoise(I, sigma);
 	I2b = convertIm(I2);
 	responses = Harris(I2, 0, R, offset);
 	compute_NMS_and_display(W, 2, I2b, responses, "noise", responses0, nms0, anms0,
-							w, h, R, offset, N, w);
+							w, h, R, offset, N, w, 0, 0, H);
 
 	// best Harris responses on scaled image
 	cout << "Scale" << endl;

TPs/move_data.sh

+# Make directory
+if [ ! -d "$2" ]; then
+	mkdir "$2"
+fi
+
+# Move files
+if ls $1/*.png 1> /dev/null 2>&1 ; then
+	ls $1/*.png | xargs mv -t "$2"
+fi
+if ls $1/*.jpeg 1> /dev/null 2>&1 ; then
+	ls $1/*.jpeg | xargs mv -t "$2"
+fi
+if ls $1/*_loc.txt 1> /dev/null 2>&1 ; then
+	ls $1/*_loc.txt | xargs mv -t "$2"
+fi
+if [ -f $1/noise.txt ]; then
+	mv $1/noise.txt -t $2
+fi
+if [ -f $1/original.txt ]; then
+	mv $1/original.txt -t $2
+fi
+if [ -f $1/rotation.txt ]; then
+	mv $1/rotation.txt -t $2
+fi
+if [ -f $1/scale.txt ]; then
+	mv $1/scale.txt -t $2
+fi
+if [ -f $1/viewpoint.txt ]; then
+	mv $1/viewpoint.txt -t $2
+fi
\ No newline at end of file

TPs/run_harris.sh

+./build/Harris im1.jpg
+./move_data.sh . im1
+python3 show.py im1/
+python3 show_loc.py im1/ 100
+python3 show_loc.py im1/ 300
+python3 show_loc.py im1/ 600
+
+./build/Harris im2.jpg
+./move_data.sh . im2
+python3 show.py im2/
+python3 show_loc.py im2/ 100
+python3 show_loc.py im2/ 300
+python3 show_loc.py im2/ 600
\ No newline at end of file

TPs/show.py

+import pylab as pl
+import argparse
+import math
+
+parser = argparse.ArgumentParser(description='Plot the number of detection in function of the threshold')
+parser.add_argument("folder", type=str, help='folder where the 12 files needed are stored')
+folder = parser.parse_args().folder
+
+trans = ["noise", "rotation", "scale", "viewpoint", "original"]
+x0, y0, x1, y1 = 1000, 1000, 0, 0
+
+for t in trans:
+	try:
+		f = open(folder+t+'.txt', "r")
+	except IOError:
+		print("Could not open the file given in parameter")
+		exit(0)
+	title = f.readline()
+	n = int(f.readline())
+	x, y = [0], [n-i for i in range(n+1)]
+	fact = 1 if t == "noise" else 5
+	tx0, tx1 = 0, 0 
+	ty0, ty1 = 0, 0 
+	for i in range(n):
+		x.append(float(f.readline()))
+		if x[-1] == x[-2]: continue
+		if x[-1] - x[-2] < fact*1e-3: tx0, ty1 = x[-1], n-i
+		elif tx1 == 0 and x[-1] - x[-2] > fact*2e3: tx1, ty0 = x[-1], n-i
+	x0, x1 = min(x0, tx0), max(x1, tx1)
+	y0, y1 = min(y0, ty0), max(y1, ty1)
+	pl.plot(x, y, label=title)
+
+pl.xlim(x0-(x1-x0)*0.02, x1)
+pl.ylim(y0-(y1-y0)*0.02, y1)
+pl.xlabel('Threshold')
+pl.ylabel('Number of detections')
+pl.legend()
+pl.savefig(f'{folder}detect.png')
+pl.show()
\ No newline at end of file

TPs/showHarris.py

-import pylab as pl
-import argparse
-import math
-
-parser = argparse.ArgumentParser(description='Plot the number of points in function of the threshold')
-parser.add_argument("file", type=str, help='file to plot')
-file = parser.parse_args().file
-
-try:
-	f = open(file, "r")
-except IOError:
-	print("Could not open the file given in parameter")
-	exit(0)
-
-title = f.readline()
-n = int(f.readline())
-x = [0]
-y = [n-i for i in range(n+1)]
-for _ in range(n):
-	x.append(float(f.readline()))
-pl.plot(x, y)
-pl.title(title)
-pl.show()
\ No newline at end of file

TPs/show_loc.py

+import pylab as pl
+import argparse
+import matplotlib as mpl
+
+mpl.rcParams['figure.subplot.left'] = 0.05
+mpl.rcParams['figure.subplot.right'] = 0.98
+mpl.rcParams['figure.subplot.bottom'] = 0.085
+mpl.rcParams['figure.subplot.top'] = 0.94
+
+parser = argparse.ArgumentParser(description='Plot the number of repeated detection in function of the threshold')
+parser.add_argument("folder", type=str, help='folder where the files needed are stored')
+parser.add_argument("N", type=str, help='Number of points')
+args = parser.parse_args()
+folder = args.folder
+N = args.N
+
+trans = ["noise", "rotation", "scale", "viewpoint"]
+alg = ["", "_nms", "_anms"]
+algNames = ["Harris", "NMS", "ANMS"]
+
+na, nt = len(alg), len(trans)
+fig, axs = pl.subplots(1, na, figsize=(6*na, 5))
+fig2, axs2 = pl.subplots(1, nt, figsize=(5*nt, 4.5))
+
+for i in range(na):
+	for t in range(nt):
+		try:
+			f = open(folder+trans[t]+alg[i]+'_'+N+'_loc.txt', "r")
+		except IOError:
+			print(f"Could not open a required file {folder+trans[t]+alg[i]+'_'+N+'_loc.txt'}")
+			exit(0)
+		title = f.readline()
+		n = int(f.readline())
+		x = [0]
+		y = [i for i in range(n+1)]
+		for _ in range(n):
+			x.append(float(f.readline()))
+		axs[i].plot(x, y, label=trans[t])
+		axs2[t].plot(x, y, label=algNames[i])
+
+for i in range(na):
+	axs[i].set_title(f'{algNames[i]} ({N} points)')
+	axs[i].set_xlabel('Threshold')
+	axs[i].set_ylabel('Number of repeated detections with original image')
+	axs[i].legend()
+for t in range(nt):
+	axs2[t].set_title(f'{trans[t]} ({N} points)')
+	axs2[t].set_xlabel('Threshold')
+	axs2[t].set_ylabel('Number of repeated detections with original image')
+	axs2[t].legend()
+fig.savefig(f"{folder}loc_{N}.png")
+fig2.savefig(f"{folder}locb_{N}.png")
+pl.show()
\ No newline at end of file