diff --git a/lib/libstereo/middlebury/main.cpp b/lib/libstereo/middlebury/main.cpp
index 1d03ebb146b7afcd4f6d68082b6b989425915db2..d1e0529ccbcaa1825b994a7ab06af99795a3278c 100644
--- a/lib/libstereo/middlebury/main.cpp
+++ b/lib/libstereo/middlebury/main.cpp
@@ -171,10 +171,30 @@ void show_result(std::string name, const MiddleburyData &data, const Result &res
 	cv::imshow(data.name + " (error) - " + name, err_color);
 
 	printf("\n%s: %s\n", name.c_str(), data.name.c_str());
-	for (auto res : result.results) {
-		printf("%9.2f %%    correct (err < %.2f)\n", 100.0f * res.err_bad_nonoccl, res.threshold);
-		printf("%9.2f px   RMSE\n", res.rms_bad_nonoccl);
-	}
+	int nt = result.results.size();
+
+	printf("%8s", "");
+	for (auto res : result.results) { printf("%10.2f ", res.threshold); }
+	printf("\n");
+
+	for (auto res : result.results) { printf("%13s", "-------------"); }
+	printf("\n");
+
+	printf("%-8s", "total");
+	for (auto res : result.results) { printf("%9.2f%% ", res.err_total*100.0f); }
+	printf("\n");
+
+	printf("%-8s", "bad");
+	for (auto res : result.results) { printf("%9.2f%% ", res.err_bad*100.0f); }
+	printf("\n");
+
+	printf("%-8s", "invalid");
+	for (auto res : result.results) { printf("%9.2f%% ", res.err_invalid*100.0f); }
+	printf("\n");
+
+	printf("%-8s", "RMSE");
+	for (auto res : result.results) { printf("%8.2fpx ", res.rms_good); }
+	printf("\n");
 }
 
 void show_results( const MiddleburyData &data, const std::map<std::string, Result> &results) {
diff --git a/lib/libstereo/middlebury/middlebury.cpp b/lib/libstereo/middlebury/middlebury.cpp
index e3136325f1275ba1c1de9983ebeb24b223f71ba8..e9c607c57030d76cea584c3954935ff2d8431d51 100644
--- a/lib/libstereo/middlebury/middlebury.cpp
+++ b/lib/libstereo/middlebury/middlebury.cpp
@@ -117,56 +117,71 @@ MiddEvalCalib read_calibration(const std::string &filename) {
 	return calib;
 }
 
-MiddEvalResult evaluate(const cv::Mat &disp, const cv::Mat &gt, const cv::Mat &mask, float threshold) {
-	MiddEvalResult result {0.0f, 0.0f, 0.0f};
-	result.threshold = threshold;
-
-	cv::Mat diff1(disp.size(), CV_32FC1);
-	cv::Mat diff2;
-
-	cv::absdiff(disp, gt, diff1);
-	diff1.setTo(0.0f, gt != gt); // NaN
-
-	// ...? where do inf values come from; issue in subpixel interpolation?
-	// doesn't seem to happen when interpolation is enabled
-	diff1.setTo(0.0f, gt == INFINITY);
-
-	float n, n_gt, n_bad, err2, err2_bad;
-
-	// errors incl. occluded areas
-	n = countNonZero(disp);
-	n_gt = countNonZero(gt);
-	cv::pow(diff1, 2, diff2);
-	err2 = cv::sum(diff2)[0];
+MiddEvalResult evaluate(const cv::Mat &disp, const cv::Mat &gtdisp, const cv::Mat &mask, float threshold) {
+	MiddEvalResult result {0, 0, 0, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f};
+	if (gtdisp.type() != CV_32FC1) { throw std::exception(); }
+	if (disp.type() != CV_32FC1) { throw std::exception(); }
+	if (!mask.empty() && mask.type() != CV_8UC1) { throw std::exception(); }
 
-	result.err_total = n/n_gt;
-	result.rms_total = sqrt(err2/n);
+	using std::min;
+	using std::max;
 
-	diff2.setTo(0.0f, diff1 > threshold);
-	n_bad = countNonZero(diff1 <= threshold);
-	err2_bad = cv::sum(diff2)[0];
-	result.err_bad = n_bad/n_gt;
-	result.rms_bad = sqrt(err2_bad/n);
+	bool usemask = !mask.empty();
 
-	// errors ignoring occlusions (mask)
-	diff1.setTo(0.0f, mask != 255);
-	cv::pow(diff1, 2, diff2);
+	result.threshold = threshold;
+	int &n = result.n;
+	int &bad = result.bad;
+	int good = 0;
+	int &invalid = result.invalid;
+
+	float &rms_good = result.rms_good;
+	float &rms_bad = result.rms_bad;
+
+	float serr = 0;
+
+	// evaldisp.cpp from middlebury SDK
+
+	for (int y = 0; y < gtdisp.rows; y++) {
+	for (int x = 0; x < gtdisp.cols; x++) {
+		float gt = gtdisp.at<float>(y, x, 0);
+		if (gt == INFINITY) // unknown
+		continue;
+
+		float d = disp.at<float>(y, x, 0);
+		int valid = (d != INFINITY && d != 0.0f); // added 0.0f
+		/*if (valid) {
+			float maxd = maxdisp; // max disp range
+			d = max(0, min(maxd, d)); // clip disps to max disp range
+		}*/
+		/*if (valid && rounddisp) { d = round(d); }*/
+		float err = fabs(d - gt);
+
+		if (usemask && mask.at<uchar>(y, x, 0) != 255) {} // don't evaluate pixel
+		else {
+			n++;
+			if (valid) {
+				serr += err;
+				if (err > threshold) {
+					bad++;
+					rms_bad += err*err;
+				}
+				else {
+					rms_good += err*err;
+					good++;
+				}
+			} else {// invalid (i.e. hole in sparse disp map)
+				invalid++;
+			}
+		}
+	}
+	}
 
-	cv::Mat tmp;
-	disp.copyTo(tmp);
-	tmp.setTo(0.0f, mask != 255);
-	n = countNonZero(tmp);
-	n_gt = countNonZero(mask == 255);
-	err2 = cv::sum(diff2)[0];
-
-	result.err_nonoccl = n/n_gt;
-	result.rms_nonoccl = sqrt(err2/n);
-
-	diff2.setTo(0.0f, diff1 > threshold);
-	n_bad = countNonZero(diff1 <= threshold) - countNonZero(mask != 255);
-	err2_bad = cv::sum(diff2)[0];
-	result.err_bad_nonoccl = n_bad/n_gt;
-	result.rms_bad_nonoccl = sqrt(err2_bad/n);
+	rms_bad = sqrtf(rms_bad/float(bad));
+	rms_good = sqrtf(rms_good/float(good));
+	result.err_bad = float(bad)/float(n);
+	result.err_invalid = float(invalid)/float(n);
+	result.err_total = float(bad+invalid)/float(n);
+	float avgErr = serr / float(n - invalid);
 	return result;
 }
 
diff --git a/lib/libstereo/middlebury/middlebury.hpp b/lib/libstereo/middlebury/middlebury.hpp
index 87d9af5d307bb10782a8723cdb728841982c0ab7..3355860eec91bc4f6f25bf0ae9e6d46156a983da 100644
--- a/lib/libstereo/middlebury/middlebury.hpp
+++ b/lib/libstereo/middlebury/middlebury.hpp
@@ -3,14 +3,14 @@
 #include <opencv2/core/mat.hpp>
 
 struct MiddEvalResult {
-	float err_total;   // all pixels
-	float rms_total;   //
-	float err_nonoccl;   // non-masked pixels
-	float rms_nonoccl;   //
-	float err_bad;     // within threshold disparity from correct value
-	float rms_bad;     // RMS for pixels within threshold disparity from correct value
-	float err_bad_nonoccl;
-	float rms_bad_nonoccl;
+	int n;
+	int bad;
+	int invalid;
+	float err_bad;
+	float err_invalid;
+	float err_total;
+	float rms_bad;
+	float rms_good;
 	float threshold;
 };