slow_flow.cpp 36 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068
  1. /*
  2. * slow_flow.cpp
  3. *
  4. * Created on: Mar 7, 2016
  5. * Author: Janai
  6. */
  7. #include "configuration.h"
  8. #include <fstream>
  9. #include <stdlib.h>
  10. #include <string.h>
  11. #include <string>
  12. #include <cmath>
  13. #include <omp.h>
  14. #include <stdio.h>
  15. #include <boost/filesystem.hpp>
  16. #include <boost/regex.hpp>
  17. #include <opencv2/core.hpp>
  18. #include <opencv2/highgui.hpp>
  19. #include <opencv2/imgproc.hpp>
  20. #include "epic_flow_extended/image.h"
  21. #include "epic_flow_extended/io.h"
  22. #include "epic_flow_extended/epic.h"
  23. #include "epic_flow_extended/variational_mt.h"
  24. #include "utils/utils.h"
  25. #include "utils/parameter_list.h"
  26. // include Hamilton-Adams demosaicing
  27. extern "C"
  28. {
  29. #ifdef DMGUNTURK
  30. #include DMGUNTURK_PATH(/dmha.h)
  31. #endif
  32. }
  33. // include flowcode (middlebury devkit)
  34. #include MIDDLEBURY_PATH(/colorcode.h)
  35. #include MIDDLEBURY_PATH(/flowIO.h)
  36. // include TRWS
  37. #include TRWS_PATH(/MRFEnergy.h)
  38. void HADemosaicing(float *Output, const float *Input, int Width, int Height, int RedX, int RedY) {
  39. #ifdef DMGUNTURK
  40. HamiltonAdamsDemosaic(Output, Input, Width, Height, RedX, RedY); // Hamilton-Adams implemented by Pascal Getreuer
  41. #endif
  42. }
  43. using namespace std;
  44. using namespace cv;
  45. enum COMPARISON {GROUNDTRUTH = 0, WARPING = 1};
  46. /* show usage information */
  47. void usage(){
  48. printf("usage:\n");
  49. printf(" ./slow_flow [cfg] -overwrite -resume -deep_settings [settings] -threads -fr [select one specific adaptive frame rate] -jet [select one specific high speed flow]\n");
  50. printf("\n");
  51. }
  52. void setDefault(ParameterList& params) {
  53. // general
  54. params.insert("verbose", "0", true);
  55. params.insert("threads", "1", true);
  56. params.insert("16bit", "1", true); // set to 1 if input images 16 bit
  57. params.insert("raw", "1", true); // set to 1 if raw images
  58. params.insert("raw_weight", "1", true); // weight for raw pixel
  59. params.insert("raw_demosaicing", "1", true); // demosaicing method 0: bilinear interp, 1: hamilton adams, 2: opencv
  60. params.insert("raw_red_loc", "1,0", true); // location of first red pixel (x,y)
  61. params.insert("Jets", "1", true); // number of high speed flow estimates
  62. params.insert("adaptive", "1", true); // choose two frame rates according to 0.99 quantile
  63. params.insert("max_fps", "200", true); // frame rate of input sequence
  64. params.insert("ref_fps", "20", true); // final frame rate
  65. params.insert("scale", "1.0f", true); // scaling factor for input images
  66. params.insert("sigma", "0.0f", true); // presmoothing
  67. params.insert("deep_matching", "1", true); // set to 1 to use deep matching
  68. params.insert("dm_scale", "1.0f", true); // scaling factor for deep matching
  69. params.insert("slow_flow_method", "symmetric", true); // symmetric: symmetric window, forward: forward window
  70. params.insert("slow_flow_S", "2", true); // number of frames in window
  71. // energy function
  72. params.insert("slow_flow_dataterm", "1", true); // 0: unnormalized, 1: normalized
  73. params.insert("slow_flow_smoothing", "1", true); // 0: \phi(u_dx) + \phi(u_dy), 1: \phi(u_dx + u_symdy) + \phi(u_dy + u_symdx), 2: \phi(u_dx + u_dy)
  74. params.insert("slow_flow_alpha", "4.0f", true); // flow smoothness weight
  75. params.insert("slow_flow_gamma", "6.0f", true); // gradient constancy assumption
  76. params.insert("slow_flow_delta", "1.0f", true); // color constancy assumption weight
  77. params.insert("slow_flow_rho_0", "1", true); // weight for successive data terms
  78. params.insert("slow_flow_rho_1", "1", true); // weight for successive data terms
  79. params.insert("slow_flow_omega_0", "0", true); // weight for reference data terms
  80. params.insert("slow_flow_omega_1", "2", true); // weight for reference data terms
  81. // image pyramid
  82. params.insert("slow_flow_layers", "1", true); // number of pyramid layers
  83. params.insert("slow_flow_p_scale", "0.9f", true); // scaling factor for pyramid
  84. // optimization
  85. params.insert("slow_flow_niter_alter", "10", true); // number of alternations
  86. params.insert("slow_flow_niter_graphc", "10", true); // number of iterations for graph cut expansion algorithm
  87. params.insert("slow_flow_niter_outer", "10", true); // number of outer fixed point iterations
  88. params.insert("slow_flow_thres_outer", "1e-5", true); // threshold for du and dv to stop the optimization
  89. params.insert("slow_flow_niter_inner", "1", true); // number of inner fixed point iterations
  90. params.insert("slow_flow_thres_inner", "1e-5", true); // threshold for du and dv to stop the optimization
  91. params.insert("slow_flow_niter_solver", "30", true); // number of solver iterations
  92. params.insert("slow_flow_sor_omega", "1.9f", true); // omega parameter of sor method
  93. // occlusion reasoning
  94. params.insert("slow_flow_occlusion_reasoning", "1", true); // set to 1 to enable occlusion reasoning
  95. params.insert("slow_flow_occlusion_penalty", "0.1", true); // preference of backwards occlusion (using the forward data terms)
  96. params.insert("slow_flow_occlusion_alpha", "0.1", true); // occlusion smoothness weight
  97. params.insert("slow_flow_output_occlusions", "1", true); // set to 1 to output occlusion estimate
  98. // regularization
  99. params.insert("slow_flow_robust_color", "1", true); // 0: quadratic, 1: modified_l1_norm, 2: lorentzian, 3: trunc mod l1, 4: Geman McClure
  100. params.insert("slow_flow_robust_color_eps", "0.001", true); // epsilon of robust function
  101. params.insert("slow_flow_robust_color_truncation", "0.5", true); // truncation threshold (trun_mod_l1)
  102. params.insert("slow_flow_robust_reg", "1", true); // 0: quadratic, 1: modified_l1_norm, 2: lorentzian, 3: trunc mod l1, 4: Geman McClure
  103. params.insert("slow_flow_robust_reg_eps", "0.001", true); // epsilon of robust function
  104. params.insert("slow_flow_robust_reg_truncation", "0.5", true); // truncation threshold (trun_mod_l1)
  105. }
  106. inline bool insideImg(double x, double y, int width, int height) {
  107. return (y >= 0 && y < height && x >= 0 && x < width);
  108. }
  109. int main(int argc, char **argv){
  110. if( argc < 2){
  111. if(argc>1) fprintf(stderr,"Error, not enough arguments\n");
  112. usage();
  113. exit(1);
  114. }
  115. int n_thread = 1;
  116. // read optional arguments
  117. string sequence_path = "", output = "", format = "";
  118. uint32_t start = 0;
  119. ParameterList params;
  120. setDefault(params); // set default parameters
  121. // read in ParameterList
  122. if( argc > 1 && argv[1][0] != '-' && access( argv[1], F_OK ) != -1) {
  123. params.read(argv[1]);
  124. } else {
  125. cerr << "Couldn't find " << argv[1] << "!" << endl;
  126. return -1;
  127. }
  128. #define isarg(key) !strcmp(a,key)
  129. bool overwrite_output = false;
  130. bool resume_frame = false;
  131. string input_path = "";
  132. string output_path = "";
  133. string deep_settings = "";
  134. double max_flow_scale = 3.0;
  135. int selected_jet = -1;
  136. int selected_fr = -1;
  137. int current_arg = 1;
  138. while(current_arg < argc ){
  139. const char* a = argv[current_arg++];
  140. if(a[0] != '-') {
  141. continue;
  142. }
  143. if( isarg("-h") || isarg("-help") )
  144. usage();
  145. else if( isarg("-overwrite") )
  146. overwrite_output = true;
  147. else if( isarg("-resume") )
  148. resume_frame = true;
  149. else if( isarg("-deep_settings") )
  150. deep_settings = string(argv[current_arg++]);
  151. else if( isarg("-threads") )
  152. params.insert("threads", argv[current_arg++], true);
  153. else if( isarg("-fr") )
  154. selected_fr = atoi(argv[current_arg++]);
  155. else if( isarg("-jet") ) {
  156. selected_jet = atoi(argv[current_arg++]);
  157. resume_frame = true;
  158. }else{
  159. fprintf(stderr, "unknown argument %s", a);
  160. usage();
  161. }
  162. }
  163. bool enable_dm = params.parameter<bool>("deep_matching");
  164. // to restrict deep call -deep_settings '-ngh_rad <max_flow>'
  165. float max_flow = 50;
  166. if(params.exists("max_flow"))
  167. max_flow = max(5.0f, params.parameter<float>("max_flow"));
  168. float scale = params.parameter<float>("scale","1.0");
  169. n_thread = params.parameter<int>("threads");
  170. start = params.sequence_start;
  171. // decompose sequence in batches
  172. int steps = params.parameter<int>("slow_flow_S") - 1;
  173. int ref = steps;
  174. vector<int> compression_params;
  175. compression_params.push_back(CV_IMWRITE_PNG_COMPRESSION);
  176. compression_params.push_back(0);
  177. compression_params.push_back(CV_IMWRITE_JPEG_QUALITY);
  178. compression_params.push_back(100);
  179. int max_fps = params.parameter<int>("max_fps", "1"); // sequence fps
  180. int jet_fps = max_fps;
  181. if(params.exists("jet_fps")) jet_fps = params.parameter<int>("jet_fps"); // jet estimation fps
  182. int skip = (1.0f * max_fps) / jet_fps; // number of frames to skip for jet fps
  183. cout << skip << endl;
  184. // split sequence path into path an format
  185. bool sintel = params.parameter<bool>("sintel", "0"); // specific file names (we would like to be able to distinguish frame number from 24 fps and 1008 fps)
  186. bool subframes = params.parameter<bool>("subframes", "0"); // are subframes specified
  187. int start_format = (params.file.find_last_of('/') + 1);
  188. int end_format = params.file.length() - start_format;
  189. sequence_path = params.file.substr(0,start_format);
  190. format = params.file.substr(start_format, end_format);
  191. if(sequence_path[sequence_path.length() - 1] != '/') sequence_path = sequence_path + "/";
  192. params.file = sequence_path;
  193. params.insert("format", format, true);
  194. if(sequence_path.empty() || params.output.empty())
  195. return -1;
  196. int len_format = format.find_last_of('.');
  197. string format_flow = format.substr(0,len_format);
  198. if(sintel && !subframes) {
  199. start = start * 1000;
  200. for(uint32_t i = 0; i < params.sequence_start_list.size(); i++) {
  201. params.sequence_start_list[i] = params.sequence_start_list[i] * 1000; //
  202. }
  203. }
  204. params.sequence_start = start;
  205. // MAKE SURE FOLDER IS NOT OVERWRITTEN
  206. if(!resume_frame && !overwrite_output) {
  207. string newPath = params.output;
  208. if(newPath[newPath.length() - 1] == '/') newPath.erase(newPath.length() - 1);
  209. int num = 1;
  210. while(boost::filesystem::exists(newPath)) {
  211. cerr << newPath << " already exists!" << endl;
  212. stringstream tmp;
  213. tmp << newPath << "_" << num++;
  214. newPath = tmp.str();
  215. }
  216. params.output = newPath;
  217. }
  218. if(params.output[params.output.length() - 1] != '/') params.output = params.output + "/";
  219. epic_params_t epic_params;
  220. epic_params_default(&epic_params);
  221. epic_params.pref_nn= 25;
  222. epic_params.nn= 160;
  223. epic_params.coef_kernel = 1.1f;
  224. /*
  225. * ################## read in quantil
  226. */
  227. bool adaptive = false;
  228. double quantil = 1.0;
  229. double hfr_quantil = 2.0;
  230. int hfr_rate = 1;
  231. int lfr_rate = 4;
  232. string adfr = SOURCE_PATH + "/adaptiveFR.dat";
  233. if (access(adfr.c_str(), F_OK) != -1) {
  234. char line[500];
  235. fstream adativef;
  236. adativef.open(adfr.c_str());
  237. while(adativef.getline(line, 500)) {
  238. char* split = strtok(line, "\n"); // get rid of '\n' at the end
  239. split = strtok(split, "\t");
  240. char* val = strtok(NULL, "\t");
  241. if(strcmp(split, "opt_hfr_quantil") == 0) {
  242. hfr_quantil = atof(val);
  243. }
  244. if(strcmp(split, "opt_lfr_rate") == 0) {
  245. lfr_rate = atof(val);
  246. }
  247. }
  248. adaptive = params.parameter<bool>("adaptive", "0");
  249. adativef.close();
  250. }
  251. double orig_max_flow = 0;
  252. string qfstr = (sequence_path+"/quantil.dat");
  253. if (!params.exists("max_flow") && access(qfstr.c_str(), F_OK) != -1) {
  254. char line[500];
  255. fstream quantilf;
  256. quantilf.open(qfstr.c_str());
  257. quantilf.getline(line, 500);
  258. quantil = atof(line);
  259. // use maximum or quantil
  260. if(quantilf.getline(line, 500))
  261. orig_max_flow = max_flow_scale * atof(line);
  262. else
  263. orig_max_flow = max_flow_scale * quantil;
  264. // compute frame rate
  265. if(adaptive) {
  266. int keyframes = params.parameter<float>("max_fps") / params.parameter<float>("ref_fps");
  267. if(keyframes == 0) {
  268. // exact rates
  269. hfr_rate = hfr_quantil / quantil;
  270. hfr_rate = max(1.0, round(hfr_rate)); // rounded and minimum 1
  271. lfr_rate = hfr_rate * lfr_rate; // too small quantils will have the same hfr and lfr this way!
  272. lfr_rate = hfr_rate * lfr_rate;
  273. // make sure we have the same keyframes
  274. double m = round(lfr_rate / hfr_rate);
  275. lfr_rate = hfr_rate * m;
  276. cout << hfr_rate << " " << lfr_rate << endl;
  277. } else {
  278. // with keyframes
  279. hfr_rate = max(1.0, round(hfr_quantil / quantil)); // rounded and minimum 1
  280. while(hfr_rate < keyframes && keyframes % (hfr_rate * steps) != 0)
  281. hfr_rate++;
  282. cout << "hfr_rate " << hfr_rate << endl;
  283. lfr_rate = min(keyframes, hfr_rate * lfr_rate);
  284. while((lfr_rate * steps < keyframes && (keyframes % (lfr_rate * steps) != 0 || (keyframes % (lfr_rate * steps) == 0 && (lfr_rate * steps) % (hfr_rate * steps) != 0))) ||
  285. (lfr_rate * steps >= keyframes && (lfr_rate * steps) % (hfr_rate * steps) != 0))
  286. lfr_rate++;
  287. lfr_rate = min(keyframes / steps, lfr_rate);
  288. cout << "lfr_rate " << lfr_rate << endl;
  289. }
  290. } else {
  291. // set maximum flow according to quantil
  292. max_flow = max(5.0, orig_max_flow * scale * ref * skip); // twice to make sure its big enough and at least 5 pixel
  293. }
  294. } else
  295. adaptive = false;
  296. int start_fr = 0;
  297. int end_fr = (adaptive + 1);
  298. if(selected_fr >= 0) {
  299. start_fr = selected_fr;
  300. end_fr = selected_fr + 1;
  301. }
  302. string orig_deep_settings = deep_settings;
  303. for(int adFR = start_fr; adFR < end_fr; adFR++) {
  304. ParameterList adaptCfg(params);
  305. deep_settings = orig_deep_settings;
  306. if(adaptive) {
  307. stringstream jfps;
  308. if(adFR == 0) {
  309. adaptCfg.output += "high_fr/";
  310. // set frame rate
  311. jfps << max_fps / hfr_rate;
  312. adaptCfg.insert("jet_fps", jfps.str(), true); // jet estimation fps
  313. skip = hfr_rate;
  314. // compute max flow
  315. max_flow = max(5.0, orig_max_flow * scale * ref * hfr_rate);
  316. } else {
  317. adaptCfg.output += "low_fr/";
  318. // set frame rate
  319. jfps << max_fps / lfr_rate;
  320. adaptCfg.insert("jet_fps", jfps.str(), true); // jet estimation fps
  321. skip = lfr_rate;
  322. // compute max flow
  323. max_flow = max(5.0, orig_max_flow * scale * ref * lfr_rate);
  324. }
  325. }
  326. // SMALLER RESOLUTION FOR DEEP MATCHING
  327. double dm_scale = params.parameter<float>("dm_scale","1.0");
  328. if(enable_dm && max_flow > 150) {
  329. dm_scale = 0.5 * dm_scale;
  330. max_flow = max(5.0, 0.5 * max_flow);
  331. }
  332. /*
  333. * frames = 1 for reference frame
  334. * + steps window before first frame
  335. * + (Jets - 1) * steps window for each jet
  336. * + steps window after last frame
  337. * + steps additional window for backward flow
  338. */
  339. int frames = 1 + (adaptCfg.Jets + 2) * steps;
  340. // only read necessary frames
  341. uint32_t start_f = 0;
  342. uint32_t end_f = frames;
  343. uint32_t start_j = 0;
  344. uint32_t end_j = adaptCfg.Jets;
  345. if(resume_frame && selected_jet >= 0) {
  346. start_f = selected_jet * steps;
  347. end_f = min(frames, 1 + (selected_jet + 3) * steps);
  348. start_j = selected_jet;
  349. end_j = min((int) adaptCfg.Jets, selected_jet + 1);
  350. }
  351. if(start_f > end_f) continue;
  352. // create results folder
  353. boost::filesystem::create_directories(adaptCfg.output); // ParameterList result folder
  354. boost::filesystem::create_directories(adaptCfg.output+"/sequence/"); // ParameterList result folder
  355. if(!adaptCfg.file_gt.empty()) boost::filesystem::create_directories(adaptCfg.output+"/gt/"); // ParameterList result folder
  356. int width = 0;
  357. int height = 0;
  358. /*
  359. * ################### read in image sequence ###################
  360. */
  361. vector<int> red_loc = adaptCfg.splitParameter<int>("raw_red_loc","0,0");
  362. char** img_files = new char*[frames];
  363. color_image_t **un_seq = new color_image_t*[frames];
  364. color_image_t **un_seq_back = new color_image_t*[frames];
  365. color_image_t **seq = new color_image_t*[frames];
  366. color_image_t **seq_back = new color_image_t*[frames];
  367. for(uint32_t f = start_f; f < end_f; f++) {
  368. char img_file[200];
  369. if(!sintel) {
  370. sprintf(img_file, (sequence_path+format).c_str(), start - ref * skip + f * skip);
  371. } else {
  372. int sintel_frame = start / 1000;
  373. int hfr_frame = f * skip - ref * skip + (start % 1000);
  374. while(hfr_frame < 0) {
  375. sintel_frame--;
  376. hfr_frame = 42 + hfr_frame;
  377. }
  378. while(hfr_frame > 41) {
  379. sintel_frame++;
  380. hfr_frame = hfr_frame - 42;
  381. }
  382. sprintf(img_file, (sequence_path+format).c_str(), sintel_frame, hfr_frame);
  383. }
  384. cout << "Reading " << img_file << "..." << endl;
  385. Mat img = imread(string(img_file), CV_LOAD_IMAGE_UNCHANGED); // load images
  386. float norm = 1;
  387. if(img.type() == 2 || img.type() == 18)
  388. norm = 1.0f/255; // for 16 bit images
  389. // convert to floating point
  390. img.convertTo(img, CV_32FC(img.channels()));
  391. /*
  392. * DEMOSAICING
  393. */
  394. if(adaptCfg.exists("raw") && adaptCfg.parameter<bool>("raw")) {
  395. Mat tmp = img.clone();
  396. color_image_t* tmp_in = color_image_new(img.cols, img.rows);
  397. color_image_t* tmp_out = color_image_new(img.cols, img.rows);
  398. switch(adaptCfg.parameter<int>("raw_demosaicing", "0")) {
  399. case 0: // use bilinear demosaicing
  400. img = Mat::zeros(tmp.rows, tmp.cols, CV_32FC3);
  401. bayer2rgbGR(tmp, img, red_loc[0], red_loc[1]); // red green
  402. break;
  403. case 1: // use hamilton adams demosaicing
  404. mat2colorImg<float>(img, tmp_in);
  405. HADemosaicing(tmp_out->c1, tmp_in->c1, tmp_in->width, tmp_in->height, red_loc[0], red_loc[1]); // Hamilton-Adams implemented by Pascal Getreuer
  406. img = Mat::zeros(img.rows, img.cols, CV_32FC3);
  407. colorImg2colorMat<Vec3f>(tmp_out, img);
  408. break;
  409. case 2: // use opencv demosaicing
  410. tmp.convertTo(tmp, CV_8UC1);
  411. img = Mat::zeros(tmp.rows, tmp.cols, CV_8UC3);
  412. int code = CV_BayerBG2RGB;
  413. if(red_loc[1] == 0) // y
  414. if(red_loc[0] == 0) // x
  415. code = CV_BayerBG2RGB;
  416. else
  417. code = CV_BayerGB2RGB;
  418. else
  419. if(red_loc[0] == 0) // x
  420. code = CV_BayerGR2RGB;
  421. else
  422. code = CV_BayerRG2RGB;
  423. cv::cvtColor(tmp, img, code); // components from second row, second column !!!!!!!!!!!!!!!!!
  424. img.convertTo(img, CV_32FC(img.channels()));
  425. break;
  426. }
  427. color_image_delete(tmp_in);
  428. color_image_delete(tmp_out);
  429. } else {
  430. // covert to RGB
  431. cv::cvtColor(img, img, CV_BGR2RGB);
  432. }
  433. if(!adaptCfg.exists("raw") || adaptCfg.parameter<float>("raw_weight", "1.0") == 1.0) {
  434. // use only a part of the images
  435. if(adaptCfg.extent.x > 0 || adaptCfg.extent.y > 0) {
  436. img = img.rowRange(Range(adaptCfg.center.y - adaptCfg.extent.y/2,adaptCfg.center.y + adaptCfg.extent.y/2));
  437. img = img.colRange(Range(adaptCfg.center.x - adaptCfg.extent.x/2,adaptCfg.center.x + adaptCfg.extent.x/2));
  438. }
  439. // rescale image with gaussian blur to avoid anti-aliasing
  440. if(scale != 1) {
  441. GaussianBlur(img, img, Size(),1/sqrt(2*scale),1/sqrt(2*scale),BORDER_REPLICATE);
  442. resize(img, img, Size(0,0), scale, scale, INTER_LINEAR);
  443. }
  444. }
  445. // print to file
  446. img_files[f] = new char[500];
  447. sprintf(img_files[f], (adaptCfg.output+"sequence/frame_%i.png").c_str(), start - ref * skip + f * skip);
  448. Mat output_img;
  449. if(adaptCfg.verbosity(WRITE_FILES)) {
  450. if(adaptCfg.parameter<bool>("16bit", "0")) {
  451. img.convertTo(output_img, CV_16UC(img.channels()));
  452. } else {
  453. img.convertTo(output_img, CV_8UC(img.channels()), norm);
  454. }
  455. cv::cvtColor(output_img, output_img, CV_RGB2BGR); // OpenCV uses BGR
  456. imwrite(img_files[f], output_img, compression_params);
  457. }
  458. width = img.cols;
  459. height = img.rows;
  460. // copy data
  461. seq[f] = color_image_new(width, height);
  462. if(img.channels() == 1) {
  463. mat2colorImg<float>(img, seq[f]);
  464. } else
  465. colorMat2colorImg<Vec3f>(img, seq[f]);
  466. // resize and copy data for deep match
  467. if(dm_scale != 1) {
  468. GaussianBlur(img, img,Size(),1/sqrt(2*dm_scale),1/sqrt(2*dm_scale),BORDER_REPLICATE);
  469. resize(img, img, Size(0,0), dm_scale, dm_scale, INTER_LINEAR);
  470. }
  471. sprintf(img_files[f], (adaptCfg.output+"sequence/frame_epic_%i.png").c_str(), start - ref * skip + f * skip);
  472. img.convertTo(img, CV_8UC(img.channels()), norm);
  473. if(f % steps == 0) {
  474. cv::cvtColor(img, output_img, CV_RGB2BGR); // OpenCV uses BGR
  475. imwrite(img_files[f], output_img, compression_params);
  476. }
  477. un_seq[f] = color_image_new(width*dm_scale, height*dm_scale);
  478. if(img.channels() == 1) {
  479. mat2colorImg<uchar>(img, un_seq[f]);
  480. } else
  481. colorMat2colorImg<Vec3b>(img, un_seq[f]);
  482. un_seq_back[frames - 1 - f] = un_seq[f];
  483. seq_back[frames - 1 - f] = seq[f];
  484. }
  485. /*
  486. * DEMOSAICING AND CHANNEL WEIGHTING
  487. */
  488. color_image_t* channel_weights = color_image_new(seq[start_f]->width, seq[start_f]->height);
  489. fill_n(channel_weights->c1, 3*channel_weights->height*channel_weights->stride, 1.0); // set channel weights to 1
  490. if(adaptCfg.exists("raw") && adaptCfg.parameter<bool>("raw"))
  491. rawWeighting(channel_weights, red_loc[0], red_loc[1], adaptCfg.parameter<float>("raw_weight", "1.0"));
  492. /*
  493. * ################### read in ground truth ###################
  494. */
  495. Mat* gt = new Mat[adaptCfg.Jets];
  496. for(u_int32_t j = start_j; j < end_j; j++) {
  497. char path[200];
  498. if(!sintel)
  499. sprintf(path, adaptCfg.file_gt.c_str(), start + j*steps);
  500. else {
  501. int sintel_frame = start / 1000;
  502. int hfr_frame = j*steps + (start % 1000);
  503. while(hfr_frame < 0) {
  504. sintel_frame--;
  505. hfr_frame = 42 + hfr_frame;
  506. }
  507. while(hfr_frame > 41) {
  508. sintel_frame++;
  509. hfr_frame = hfr_frame - 42;
  510. }
  511. sprintf(path, adaptCfg.file_gt.c_str(), sintel_frame, hfr_frame);
  512. }
  513. cout << path << endl;
  514. if(access( path, F_OK ) != -1) {
  515. gt[j] = readGTMiddlebury(string(path));
  516. // crop images
  517. if(adaptCfg.center.x > 0) {
  518. // NOTE: ROWRANGE IS INDUCING ERRORS IN ACCUMULATION!!!!!!
  519. gt[j] = crop(gt[j], adaptCfg.center, adaptCfg.extent);
  520. }
  521. // rescale image
  522. resize(gt[j], gt[j], Size(0, 0), scale, scale, INTER_NEAREST); // LINEAR PROBLEMATIC AT MOTION DISCONTINOUTIES
  523. gt[j] *= scale;
  524. Mat img = flowColorImg(gt[j], adaptCfg.verbosity(VER_CMD));
  525. if(! img.data) { // Check for invalid input
  526. cout << "No gt flow for frame " << endl ;
  527. continue;
  528. }
  529. // write flow and flow image to file
  530. if(!gt[j].empty() && !adaptCfg.output.empty()) {
  531. char gtF[200];
  532. sprintf(gtF, "%s/gt/flow_%05i.png", adaptCfg.output.c_str(), adaptCfg.sequence_start + j*steps);
  533. imwrite(gtF, img);
  534. }
  535. char gtF[200];
  536. sprintf(gtF, "%s/gt/flow_%05i.flo", adaptCfg.output.c_str(), adaptCfg.sequence_start + j*steps);
  537. writeFlowMiddlebury(gt[j], gtF);
  538. // DEBUG: show groundtruth flow
  539. if(adaptCfg.verbosity(VER_IN_GT)) {
  540. stringstream title;
  541. title << "GT flow of frame " << adaptCfg.sequence_start + j*steps;
  542. namedWindow( title.str(), WINDOW_FREERATIO ); // Create a window for display.
  543. imshow( title.str(), img ); // Show our image inside it.
  544. waitKey(0);
  545. }
  546. }
  547. }
  548. // normalize intensities
  549. normalize(&seq[start_f], end_f - start_f, adaptCfg);
  550. boost::filesystem::create_directories(adaptCfg.output+"tmp/"); // temporary folder for edges and deep_matches
  551. if(adaptCfg.parameter<bool>("slow_flow_occlusion_reasoning", "0"))
  552. boost::filesystem::create_directories(adaptCfg.output+"occlusion/"); // folder for occlusion output
  553. /*
  554. * write infos to file
  555. */
  556. adaptCfg.print();
  557. ofstream infos;
  558. infos.open((adaptCfg.output + "config.cfg").c_str());
  559. infos << "# SlowFlow variational estimation\n";
  560. infos << adaptCfg;
  561. infos.close();
  562. // write stats
  563. stringstream results;
  564. results << "frame\ttime\n\n";
  565. int avg_time = 0;
  566. int counter = 0;
  567. if(enable_dm && max_flow < 300) {
  568. stringstream mfstr;
  569. mfstr << ceil(max_flow);
  570. deep_settings = " -ngh_rad " + mfstr.str();
  571. cout << deep_settings << endl;
  572. cout << "Max flow: " << max_flow << endl;
  573. } else
  574. deep_settings = "";
  575. #pragma omp parallel for num_threads(n_thread) shared(seq, seq_back, gt)
  576. for(uint32_t j = start_j; j < end_j; j++) {
  577. ParameterList thread_params(adaptCfg);
  578. int f = j*steps;
  579. char curr_f [33];
  580. sprintf(curr_f, "%d_forward", start + j*steps) ;
  581. thread_params.insert("current_frame", curr_f, true);
  582. /*
  583. * indices 0 1 2 3 4 5 6 7 8 9 F-4 F-3 F-2 F-1 F F+1 F+2 F+3 F+4 F+5
  584. * forward (S = 3) : 0 1 2 (3->4->5->6) 7 8 9 ...
  585. * backward (S = 3): ... 9 8 7 (6-> 5-> 4-> 3) 2 1 0
  586. */
  587. color_image_t **un_im = &un_seq[f];
  588. color_image_t **un_im_back = &un_seq_back[frames - 1 - f - 3 * steps];
  589. color_image_t **im = &seq[f];
  590. color_image_t **im_back = &seq_back[frames - 1 - f - 3 * steps]; // -1 because zero index and -steps because different reference frame
  591. // prepare variables
  592. image_t *wx, *wy;
  593. char edges_f[1000], edges_b[1000], edges_cmd[1000], match_f[1000], match_b[1000], match_cmd[1000];
  594. time_t t_start, t_stop;
  595. time_t pp_start, pp_stop;
  596. time_t ep_start, ep_stop;
  597. int t_preprocessing = 0;
  598. /*
  599. * Compute EDGES AND MATCHES!!!!!
  600. */
  601. sprintf(edges_f, (adaptCfg.output+"tmp/edges_%i.dat").c_str(), adaptCfg.sequence_start + f * skip);
  602. sprintf(edges_b, (adaptCfg.output+"tmp/edges_%i.dat").c_str(), adaptCfg.sequence_start + f * skip + ref * skip);
  603. sprintf(match_f, (adaptCfg.output+"tmp/matches_%i_%i.dat").c_str(), adaptCfg.sequence_start + f * skip, adaptCfg.sequence_start + f * skip + ref * skip);
  604. sprintf(match_b, (adaptCfg.output+"tmp/matches_%i_%i.dat").c_str(), adaptCfg.sequence_start + f * skip + ref * skip, adaptCfg.sequence_start + f * skip);
  605. if(enable_dm) {
  606. if(!resume_frame || (resume_frame && access( edges_f, F_OK ) == -1)) {
  607. cout << "Computing edges ..." << endl;
  608. sprintf(edges_cmd, "matlab -nodesktop -nojvm -r \"addpath(\'%s/matlab/\'); detect_edges(\'%s\',\'%s\'); exit\"", SOURCE_PATH.c_str(), img_files[j * steps + ref], edges_f);
  609. time(&pp_start);
  610. system(edges_cmd);
  611. time(&pp_stop);
  612. t_preprocessing += (int) difftime(pp_stop, pp_start);
  613. }
  614. if(!resume_frame || (resume_frame && access( edges_b, F_OK ) == -1)) {
  615. cout << "Computing edges ..." << endl;
  616. sprintf(edges_cmd, "matlab -nodesktop -nojvm -r \"addpath(\'%s/matlab/\'); detect_edges(\'%s\',\'%s\'); exit\"", SOURCE_PATH.c_str(), img_files[j * steps + 2*ref], edges_b);
  617. time(&pp_start);
  618. system(edges_cmd);
  619. time(&pp_stop);
  620. t_preprocessing += (int) difftime(pp_stop, pp_start);
  621. }
  622. if(!resume_frame || (resume_frame && access( match_f, F_OK ) == -1)) {
  623. cout << "Computing matches between " << adaptCfg.sequence_start + f * skip << " and " << adaptCfg.sequence_start + f * skip + ref * skip << " ..." << endl;
  624. sprintf(match_cmd, "%s/deepmatching \'%s\' \'%s\' -png_settings %s -out \'%s\'", DEEPMATCHING_PATH.c_str(), img_files[j * steps + ref], img_files[j * steps + 2*ref], deep_settings.c_str(), match_f);
  625. time(&pp_start);
  626. system(match_cmd);
  627. time(&pp_stop);
  628. t_preprocessing += (int) difftime(pp_stop, pp_start);
  629. }
  630. if(!resume_frame || (resume_frame && access( match_b, F_OK ) == -1)) {
  631. cout << "Computing matches between " << adaptCfg.sequence_start + f * skip + ref * skip << " and " << adaptCfg.sequence_start + f * skip<< " ..." << endl;
  632. sprintf(match_cmd, "%s/deepmatching \'%s\' \'%s\' -png_settings %s -out \'%s\'", DEEPMATCHING_PATH.c_str(), img_files[j * steps + 2*ref], img_files[j * steps + ref], deep_settings.c_str(), match_b);
  633. time(&pp_start);
  634. system(match_cmd);
  635. time(&pp_stop);
  636. t_preprocessing += (int) difftime(pp_stop, pp_start);
  637. }
  638. }
  639. char forward_flow_file[200];
  640. if(!sintel)
  641. sprintf(forward_flow_file, (adaptCfg.output + format_flow + ".flo").c_str(), start + f * skip);
  642. else
  643. sprintf(forward_flow_file, (adaptCfg.output + format_flow + ".flo").c_str(), start + f * skip, 0);
  644. // skip finished frames
  645. if(!resume_frame || (resume_frame && access( forward_flow_file, F_OK ) == -1)) {
  646. /*
  647. * ################### extract edges and get matches ###################
  648. */
  649. Mat floImg;
  650. // use deep matching instead of pyramid
  651. if(enable_dm) {
  652. wx = image_new(im[ref]->width*dm_scale, im[ref]->height*dm_scale);
  653. wy = image_new(im[ref]->width*dm_scale, im[ref]->height*dm_scale);
  654. image_erase(wx);
  655. image_erase(wy);
  656. // matches to target frame
  657. float_image forward_edges = read_edges(edges_f, un_im[ref]->width, un_im[ref]->height);
  658. time(&pp_start);
  659. float_image forward_matches = read_matches(match_f);
  660. time(&pp_stop);
  661. t_preprocessing += (int) difftime(pp_stop, pp_start);
  662. color_image_t *imlab = rgb_to_lab(un_im[ref]);
  663. // initilize with deep matches
  664. cout << "Epic interpolation of forward flow ..." << endl;
  665. time(&ep_start);
  666. epic(wx, wy, imlab, &forward_matches, &forward_edges, &epic_params, 1);
  667. time(&ep_stop);
  668. // rescale flow
  669. float fx = im[ref]->width / wx->width;
  670. float fy = im[ref]->height / wx->height;
  671. if(fx != 1) {
  672. Mat tmpx(wx->height, wx->width, CV_32FC1);
  673. Mat tmpy(wy->height, wy->width, CV_32FC1);
  674. img2mat<float>(wx, tmpx);
  675. img2mat<float>(wy, tmpy);
  676. resize(tmpx, tmpx, Size(im[ref]->width, im[ref]->height), 0, 0, INTER_LINEAR);
  677. resize(tmpy, tmpy, Size(im[ref]->width, im[ref]->height), 0, 0, INTER_LINEAR);
  678. image_delete(wx); image_delete(wy);
  679. wx = image_new(im[ref]->width, im[ref]->height);
  680. wy = image_new(im[ref]->width, im[ref]->height);
  681. mat2img<float>(tmpx, wx);
  682. mat2img<float>(tmpy, wy);
  683. }
  684. image_mul_scalar(wx, fx / steps); // scale flow vectors
  685. image_mul_scalar(wy, fy / steps); // scale flow vectors
  686. floImg = flowColorImg(wx, wy, adaptCfg.verbosity(VER_CMD));
  687. // write flow image to file
  688. if(adaptCfg.verbosity(WRITE_FILES)) {
  689. if(!floImg.data) { // Check for invalid input
  690. cout << "No forward flow for frame " << start + f * skip << std::endl ;
  691. } else {
  692. if(!adaptCfg.output.empty()) {
  693. stringstream flowF;
  694. flowF << adaptCfg.output << "tmp/frame_" << start + f * skip << "_INIT.png";
  695. imwrite((flowF.str()), floImg);
  696. }
  697. }
  698. }
  699. color_image_delete(imlab);
  700. free(forward_matches.pixels);
  701. free(forward_edges.pixels);
  702. } else {
  703. wx = image_new(im[ref]->width, im[ref]->height);
  704. wy = image_new(im[ref]->width, im[ref]->height);
  705. image_erase(wx);
  706. image_erase(wy);
  707. }
  708. /*
  709. * ############ forward flow ##################
  710. */
  711. cout << "Forward flow estimation ..." << endl;
  712. Variational_MT minimzer_f;
  713. minimzer_f.setChannelWeights(channel_weights);
  714. if(thread_params.verbosity(WRITE_FILES) && thread_params.parameter<bool>("slow_flow_output_occlusions","0")) {
  715. if(thread_params.parameter<bool>("slow_flow_occlusion_reasoning","0")) {
  716. stringstream occF;
  717. occF << adaptCfg.output << "tmp/frame_" << start + f * skip << "_";
  718. thread_params.insert("slow_flow_occlusions_output", occF.str(), true);
  719. }
  720. }
  721. // energy minimization
  722. time(&t_start);
  723. minimzer_f.variational(wx, wy, im, thread_params);
  724. time(&t_stop);
  725. // output the occlusions
  726. if(thread_params.parameter<bool>("slow_flow_output_occlusions","0")) {
  727. image_t* occlusions = minimzer_f.getOcclusions();
  728. Mat occ_mat(im[ref]->height, im[ref]->width, CV_32FC1);
  729. img2mat<float>(occlusions, occ_mat);
  730. occ_mat = 0.5 * (occ_mat + 1);
  731. occ_mat.convertTo(occ_mat, CV_8UC1, 255);
  732. stringstream occF;
  733. occF << adaptCfg.output << "/occlusion/frame_" << start + f * skip << ".pbm";
  734. vector<int> compression_params_occ;
  735. compression_params_occ.push_back(CV_IMWRITE_PXM_BINARY);
  736. compression_params_occ.push_back(1); // store as binary image
  737. imwrite(occF.str(), occ_mat, compression_params_occ);
  738. }
  739. // write output file
  740. image_mul_scalar(wx, steps); // scale flow
  741. image_mul_scalar(wy, steps); // scale flow
  742. writeFlowFile(forward_flow_file, wx, wy);
  743. floImg = flowColorImg(wx, wy, adaptCfg.verbosity(VER_CMD));
  744. // write flow image to file
  745. if(!floImg.data) { // Check for invalid input
  746. cout << "No forward flow for frame " << start + f * skip << std::endl ;
  747. } else {
  748. if(!adaptCfg.output.empty()) {
  749. stringstream flowF;
  750. flowF << adaptCfg.output << "frame_" << start + f * skip << ".png";
  751. imwrite((flowF.str()), floImg);
  752. }
  753. }
  754. int time = (int) difftime(t_stop,t_start) + (int) difftime(ep_stop,ep_start) + t_preprocessing;
  755. // store results
  756. #pragma omp critical (sum)
  757. {
  758. avg_time += time;
  759. counter++;
  760. // add epe for this frame
  761. results << f * skip << "\t " << time << "\n";
  762. }
  763. // clean up
  764. image_delete(wx);
  765. image_delete(wy);
  766. cout << "Forward flow from frame " << start + f << " to " << start + f * skip + steps * skip << " finished! (Computation took " << time << " s)" << endl;
  767. } else
  768. cout << "Forward flow from frame " << start + f << " to " << start + f * skip + steps * skip << " already exist!" << endl;
  769. /*
  770. * ############ backward flow ##################
  771. */
  772. thread_params = ParameterList(adaptCfg);
  773. char backward_flow_file[200];
  774. if(!sintel)
  775. sprintf(backward_flow_file, (adaptCfg.output + format_flow + "_back.flo").c_str(), start + f * skip + steps * skip);
  776. else
  777. sprintf(backward_flow_file, (adaptCfg.output + format_flow + "_back.flo").c_str(), start + f * skip + steps * skip, 0);
  778. // skip finished frames
  779. if(!resume_frame || (resume_frame && access( backward_flow_file, F_OK ) == -1)) {
  780. t_preprocessing = 0;
  781. // use deep matching instead of pyramid
  782. if(enable_dm) {
  783. wx = image_new(im[ref]->width*dm_scale, im[ref]->height*dm_scale);
  784. wy = image_new(im[ref]->width*dm_scale, im[ref]->height*dm_scale);
  785. image_erase(wx);
  786. image_erase(wy);
  787. float_image backward_edges = read_edges(edges_b, un_im_back[ref]->width, un_im_back[ref]->height);
  788. time(&pp_start);
  789. float_image backward_matches = read_matches(match_b);
  790. time(&pp_stop);
  791. t_preprocessing += (int) difftime(pp_stop, pp_start);
  792. color_image_t *imlab = rgb_to_lab(un_im_back[ref]);
  793. // initilize with deep matches
  794. cout << "Epic interpolation of backward flow ..." << endl;
  795. time(&ep_start);
  796. epic(wx, wy, imlab, &backward_matches, &backward_edges, &epic_params, 1);
  797. time(&ep_stop);
  798. // rescale flow
  799. float fx = im[ref]->width / wx->width;
  800. float fy = im[ref]->height / wx->height;
  801. if(fx != 1) {
  802. Mat tmpx(wx->height, wx->width, CV_32FC1);
  803. Mat tmpy(wy->height, wy->width, CV_32FC1);
  804. img2mat<float>(wx, tmpx);
  805. img2mat<float>(wy, tmpy);
  806. resize(tmpx, tmpx, Size(im[ref]->width, im[ref]->height), 0, 0, INTER_LINEAR);
  807. resize(tmpy, tmpy, Size(im[ref]->width, im[ref]->height), 0, 0, INTER_LINEAR);
  808. image_delete(wx); image_delete(wy);
  809. wx = image_new(im[ref]->width, im[ref]->height);
  810. wy = image_new(im[ref]->width, im[ref]->height);
  811. mat2img<float>(tmpx, wx);
  812. mat2img<float>(tmpy, wy);
  813. }
  814. image_mul_scalar(wx, fx / steps); // scale flow vectors
  815. image_mul_scalar(wy, fy / steps); // scale flow vectors
  816. color_image_delete(imlab);
  817. free(backward_matches.pixels);
  818. free(backward_edges.pixels);
  819. } else {
  820. wx = image_new(im[ref]->width, im[ref]->height);
  821. wy = image_new(im[ref]->width, im[ref]->height);
  822. image_erase(wx);
  823. image_erase(wy);
  824. }
  825. // energy minimization
  826. cout << "Backward flow estimation ..." << endl;
  827. Variational_MT minimzer_b;
  828. if(adaptCfg.exists("method") && adaptCfg.parameter("method").compare("forward") == 0)
  829. minimzer_b.one_direction = true;
  830. time(&t_start);
  831. minimzer_b.variational(wx, wy, im_back, thread_params);
  832. time(&t_stop);
  833. // write output file
  834. image_mul_scalar(wx, steps); // scale flow
  835. image_mul_scalar(wy, steps); // scale flow
  836. writeFlowFile(backward_flow_file, wx, wy);
  837. int time = (int) difftime(t_stop,t_start) + (int) difftime(ep_stop,ep_start) + t_preprocessing;
  838. // store results
  839. #pragma omp critical (sum)
  840. {
  841. avg_time += time;
  842. counter++;
  843. }
  844. // clean up
  845. image_delete(wx);
  846. image_delete(wy);
  847. cout << "Backward flow from frame " << start + f * skip << " to " << start + f * skip + steps * skip << " finished! (Computation took " << time << " s)" << endl;
  848. } else
  849. cout << "Backward flow from frame " << start + f * skip << " to " << start + f * skip + steps * skip << " already exist!" << endl;
  850. }
  851. // clean up
  852. for(uint32_t f = start_f; f < end_f; f++) {
  853. color_image_delete(seq[f]);
  854. color_image_delete(un_seq[f]);
  855. delete[] img_files[f];
  856. // seq_back[f] only pointer!
  857. }
  858. color_image_delete(channel_weights);
  859. delete[] img_files;
  860. delete[] un_seq;
  861. delete[] un_seq_back;
  862. delete[] seq;
  863. delete[] seq_back;
  864. delete[] gt;
  865. }
  866. cout << "Done!" << endl;
  867. return 0;
  868. }