#include <ftl/middlebury.hpp> #include <glog/logging.h> #include <ftl/disparity.hpp> #include <string> #include <algorithm> using cv::Mat; using cv::Size; using std::string; using std::min; using std::max; static void skip_comment(FILE *fp) { // skip comment lines in the headers of pnm files char c; while ((c=getc(fp)) == '#') while (getc(fp) != '\n') ; ungetc(c, fp); } static void skip_space(FILE *fp) { // skip white space in the headers or pnm files char c; do { c = getc(fp); } while (c == '\n' || c == ' ' || c == '\t' || c == '\r'); ungetc(c, fp); } static void read_header(FILE *fp, const char *imtype, char c1, char c2, int *width, int *height, int *nbands, int thirdArg) { // read the header of a pnmfile and initialize width and height char c; if (getc(fp) != c1 || getc(fp) != c2) LOG(FATAL) << "ReadFilePGM: wrong magic code for " << imtype << " file"; skip_space(fp); skip_comment(fp); skip_space(fp); fscanf(fp, "%d", width); skip_space(fp); fscanf(fp, "%d", height); if (thirdArg) { skip_space(fp); fscanf(fp, "%d", nbands); } // skip SINGLE newline character after reading image height (or third arg) c = getc(fp); if (c == '\r') // <cr> in some files before newline c = getc(fp); if (c != '\n') { if (c == ' ' || c == '\t' || c == '\r') LOG(FATAL) << "newline expected in file after image height"; else LOG(FATAL) << "whitespace expected in file after image height"; } } // check whether machine is little endian static int littleendian() { int intval = 1; uchar *uval = (uchar *)&intval; return uval[0] == 1; } // 1-band PFM image, see http://netpbm.sourceforge.net/doc/pfm.html // 3-band not yet supported void ftl::middlebury::readFilePFM(Mat &img, const string &filename) { // Open the file and read the header FILE *fp = fopen(filename.c_str(), "rb"); if (fp == 0) LOG(FATAL) << "ReadFilePFM: could not open \"" << filename << "\""; int width, height, nBands; read_header(fp, "PFM", 'P', 'f', &width, &height, &nBands, 0); skip_space(fp); float scalef; fscanf(fp, "%f", &scalef); // scale factor (if negative, little endian) // skip SINGLE newline character after reading third arg char c = getc(fp); if (c == '\r') // <cr> in some files before newline c = getc(fp); if (c != '\n') { if (c == ' ' || c == '\t' || c == '\r') LOG(FATAL) << "newline expected in file after scale factor"; else LOG(FATAL) << "whitespace expected in file after scale factor"; } // Allocate the image if necessary img = Mat(height, width, CV_32FC1); // Set the image shape //Size sh = img.size(); int littleEndianFile = (scalef < 0); int littleEndianMachine = littleendian(); int needSwap = (littleEndianFile != littleEndianMachine); //printf("endian file = %d, endian machine = %d, need swap = %d\n", // littleEndianFile, littleEndianMachine, needSwap); for (int y = height-1; y >= 0; y--) { // PFM stores rows top-to-bottom!!!! int n = width; float* ptr = &img.at<float>(y, 0, 0); if ((int)fread(ptr, sizeof(float), n, fp) != n) LOG(FATAL) << "ReadFilePFM(" << filename << "): file is too short"; if (needSwap) { // if endianness doesn't agree, swap bytes uchar* ptr = (uchar *)&img.at<uchar>(y, 0, 0); int x = 0; uchar tmp = 0; while (x < n) { tmp = ptr[0]; ptr[0] = ptr[3]; ptr[3] = tmp; tmp = ptr[1]; ptr[1] = ptr[2]; ptr[2] = tmp; ptr += 4; x++; } } } if (fclose(fp)) LOG(FATAL) << "ReadFilePGM(" << filename << "): error closing file"; } // 1-band PFM image, see http://netpbm.sourceforge.net/doc/pfm.html // 3-band not yet supported void ftl::middlebury::writeFilePFM(const Mat &img, const char* filename, float scalefactor) { // Write a PFM file Size sh = img.size(); int nBands = img.channels(); if (nBands != 1) LOG(FATAL) << "WriteFilePFM(" << filename << "): can only write 1-band image as pfm for now"; // Open the file FILE *stream = fopen(filename, "wb"); if (stream == 0) LOG(FATAL) << "WriteFilePFM: could not open " << filename; // sign of scalefact indicates endianness, see pfms specs if (littleendian()) scalefactor = -scalefactor; // write the header: 3 lines: Pf, dimensions, scale factor (negative val == little endian) fprintf(stream, "Pf\n%d %d\n%f\n", sh.width, sh.height, scalefactor); int n = sh.width; // write rows -- pfm stores rows in inverse order! for (int y = sh.height-1; y >= 0; y--) { const float* ptr = &img.at<float>(0, y, 0); if ((int)fwrite(ptr, sizeof(float), n, stream) != n) LOG(FATAL) << "WriteFilePFM(" << filename << "): file is too short"; } // close file if (fclose(stream)) LOG(FATAL) << "WriteFilePFM(" << filename << "): error closing file"; } void ftl::middlebury::evaldisp(const Mat &disp, const Mat >disp, const Mat &mask, float badthresh, int maxdisp, int rounddisp) { Size sh = gtdisp.size(); Size sh2 = disp.size(); Size msh = mask.size(); int width = sh.width, height = sh.height; int width2 = sh2.width, height2 = sh2.height; int scale = width / width2; if ((!(scale == 1 || scale == 2 || scale == 4)) || (scale * width2 != width) || (scale * height2 != height)) { printf(" disp size = %4d x %4d\n", width2, height2); printf("GT disp size = %4d x %4d\n", width, height); LOG(ERROR) << "GT disp size must be exactly 1, 2, or 4 * disp size"; } int usemask = (msh.width > 0 && msh.height > 0); if (usemask && (msh != sh)) LOG(ERROR) << "mask image must have same size as GT"; int n = 0; int bad = 0; int invalid = 0; float serr = 0; for (int y = 0; y < height; y++) { for (int x = 0; x < width; x++) { float gt = gtdisp.at<float>(y, x, 0); if (gt == INFINITY) // unknown continue; float d = scale * disp.at<float>(y / scale, x / scale, 0); int valid = (!isnanf(d) && d < 256.0f); // NOTE: Is meant to be infinity in middlebury if (valid) { float maxd = scale * maxdisp; // max disp range d = max(0.0f, min(maxd, d)); // clip disps to max disp range } if (valid && rounddisp) d = round(d); float err = fabs(d - gt); if (usemask && mask.at<float>(y, x, 0) != 255) { // don't evaluate pixel } else { n++; if (valid) { serr += err; if (err > badthresh) { bad++; } } else {// invalid (i.e. hole in sparse disp map) invalid++; } } } } float badpercent = 100.0*bad/n; float invalidpercent = 100.0*invalid/n; float totalbadpercent = 100.0*(bad+invalid)/n; float avgErr = serr / (n - invalid); // CHANGED 10/14/2014 -- was: serr / n printf("mask bad%.1f invalid totbad avgErr\n", badthresh); printf("%4.1f %6.2f %6.2f %6.2f %6.2f\n", 100.0*n/(width * height), badpercent, invalidpercent, totalbadpercent, avgErr); } void ftl::middlebury::test(nlohmann::json &config) { // Load dataset images Mat l = cv::imread((string)config["middlebury"]["dataset"] + "/im0.png"); Mat r = cv::imread((string)config["middlebury"]["dataset"] + "/im1.png"); // Load ground truth Mat gt; readFilePFM(gt, (string)config["middlebury"]["dataset"] + "/disp0.pfm"); if ((float)config["middlebury"]["scale"] != 1.0f) { float scale = (float)config["middlebury"]["scale"]; //cv::resize(gt, gt, cv::Size(gt.cols * scale,gt.rows * scale), 0, 0, cv::INTER_LINEAR); cv::resize(l, l, cv::Size(l.cols * scale,l.rows * scale), 0, 0, cv::INTER_LINEAR); cv::resize(r, r, cv::Size(r.cols * scale,r.rows * scale), 0, 0, cv::INTER_LINEAR); } // Run algorithm auto disparity = ftl::Disparity::create(config["disparity"]); Mat disp; disparity->compute(l,r,disp); //disp.convertTo(disp, CV_32F); // Display results evaldisp(disp, gt, Mat(), (float)config["middlebury"]["threshold"], (int)config["disparity"]["maximum"], 0); /*if (gt.cols > 1600) { cv::resize(gt, gt, cv::Size(gt.cols * 0.25,gt.rows * 0.25), 0, 0, cv::INTER_LINEAR); }*/ if (disp.cols > 1600) { cv::resize(disp, disp, cv::Size(disp.cols * 0.25,disp.rows * 0.25), 0, 0, cv::INTER_LINEAR); } cv::resize(gt, gt, cv::Size(disp.cols,disp.rows), 0, 0, cv::INTER_LINEAR); double mindisp, mindisp_gt; double maxdisp, maxdisp_gt; Mat mask; threshold(disp,mask,255.0, 255, cv::THRESH_BINARY_INV); normalize(mask, mask, 0, 255, cv::NORM_MINMAX, CV_8U); cv::minMaxLoc(disp, &mindisp, &maxdisp, 0, 0, mask); cv::minMaxLoc(gt, &mindisp_gt, &maxdisp_gt, 0, 0); //disp = (disp < 256.0f); //disp = disp + (mindisp_gt - mindisp); disp.convertTo(disp, CV_8U, 255.0f / (maxdisp_gt*(float)config["middlebury"]["scale"])); disp = disp & mask; gt = gt / maxdisp_gt; // TODO Read from calib.txt gt.convertTo(gt, CV_8U, 255.0f); //disp = disp / maxdisp; imshow("Ground Truth", gt); imshow("Disparity", disp); imshow("Diff", gt - disp); while (cv::waitKey(10) != 27); /*cv::putText(yourImageMat, "Here is some text", cv::Point(5,5), // Coordinates cv::FONT_HERSHEY_COMPLEX_SMALL, // Font 1.0, // Scale. 2.0 = 2x bigger cv::Scalar(255,255,255), // BGR Color 1, // Line Thickness (Optional) cv::CV_AA); // Anti-alias (Optional)*/ }