1 // This is brl/bseg/boxm2/pro/processes/boxm2_roi_init_geotiff_process.cxx
2 #include <cmath>
3 #include <fstream>
4 #include <ios>
5 #include <iostream>
6 #include <algorithm>
7 #include <bprb/bprb_func_process.h>
8 //:
9 // \file
10 // \brief A process to retrieve geocam from a geotiff image and lvcs from the scene
11 // Also creates a generic camera
12 // and crops the portion of the scene for the given scene and return it
13 //
14 // \author Ozge C. Ozcanli
15 // \date May 07, 2012
16
17 #include <boxm2/boxm2_scene.h>
18
19 #include <vpgl/file_formats/vpgl_geo_camera.h>
20 #include <vpgl/algo/vpgl_camera_convert.h>
21 #include "vil/vil_image_view.h"
22 #include "vil/vil_image_resource.h"
23 #include <vil/file_formats/vil_geotiff_header.h>
24 #include <vil/file_formats/vil_tiff.h>
25 #include "vil/vil_resample_bilin.h"
26 #include "vil/vil_load.h"
27 #ifdef _MSC_VER
28 # include "vcl_msvc_warnings.h"
29 #endif
30
31 namespace boxm2_roi_init_geotiff_process_globals
32 {
33 constexpr unsigned n_inputs_ = 4;
34 constexpr unsigned n_outputs_ = 2;
35 }
36
boxm2_roi_init_geotiff_process_cons(bprb_func_process & pro)37 bool boxm2_roi_init_geotiff_process_cons(bprb_func_process& pro)
38 {
39 using namespace boxm2_roi_init_geotiff_process_globals;
40
41 //process takes 4 inputs
42 std::vector<std::string> input_types_(n_inputs_);
43 input_types_[0] = "boxm2_scene_sptr";
44 input_types_[1] = "vpgl_camera_double_sptr"; // geocam
45 input_types_[2] = "vcl_string"; // tiff image
46 input_types_[3] = "unsigned"; // resolution level, if 0 orig img res, 1 => downsample by 2^1, 2 => by 2^2
47
48 // process has 2 outputs:
49 std::vector<std::string> output_types_(n_outputs_);
50 output_types_[0] = "vpgl_camera_double_sptr"; // output generic cam
51 output_types_[1] = "vil_image_view_base_sptr";
52
53 return pro.set_input_types(input_types_) && pro.set_output_types(output_types_);
54 }
55
boxm2_roi_init_geotiff_process(bprb_func_process & pro)56 bool boxm2_roi_init_geotiff_process(bprb_func_process& pro)
57 {
58 using namespace boxm2_roi_init_geotiff_process_globals;
59
60 if ( pro.n_inputs() < n_inputs_ ) {
61 std::cout << pro.name() << ": The input number should be " << n_inputs_<< std::endl;
62 return false;
63 }
64 //get the inputs
65 boxm2_scene_sptr scene = pro.get_input<boxm2_scene_sptr>(0);
66 vpgl_camera_double_sptr cam = pro.get_input<vpgl_camera_double_sptr>(1);
67 std::string geotiff_fname = pro.get_input<std::string>(2);
68 auto level = pro.get_input<unsigned>(3);
69
70 vil_image_resource_sptr img_res = vil_load_image_resource(geotiff_fname.c_str());
71 vpgl_geo_camera* geocam = nullptr;
72 if (cam) {
73 std::cout << "Using the provided loaded camera.\n";
74 geocam = dynamic_cast<vpgl_geo_camera*> (cam.ptr());
75 } else {
76 std::cout << "Using the camera initialized from the image.\n";
77 vpgl_lvcs_sptr lvcs = new vpgl_lvcs(scene->lvcs());
78 vpgl_geo_camera::init_geo_camera(img_res, lvcs, geocam);
79 }
80
81 if (!geocam) {
82 std::cerr << "In boxm2_roi_init_geotiff_process() - the geocam could not be initialized!\n";
83 return false;
84 }
85
86 double e1,n1,e2,n2,elev,lat,lon;
87 scene->lvcs().get_origin(lat, lon, elev);
88 geocam->img_four_corners_in_utm(img_res->ni(), img_res->nj(), elev, e1, n1, e2, n2);
89 std::cout << "img res ni: " << img_res->ni() << " nj: " << img_res->nj() << std::endl;
90
91 vgl_box_3d<double> scene_bbox = scene->bounding_box();
92 std::cout.setf(std::ios::fixed, std::ios::floatfield);
93
94 // crop the image
95 vgl_box_2d<double> proj_bbox;
96 double u,v;
97 geocam->project(scene_bbox.min_x(), scene_bbox.min_y(), scene_bbox.min_z(), u, v);
98 proj_bbox.add(vgl_point_2d<double>(u,v));
99 geocam->project(scene_bbox.max_x(), scene_bbox.max_y(), scene_bbox.max_z(), u, v);
100 proj_bbox.add(vgl_point_2d<double>(u,v));
101 int min_i = int(std::max(0.0, std::floor(proj_bbox.min_x())));
102 int min_j = int(std::max(0.0, std::floor(proj_bbox.min_y())));
103 int max_i = int(std::min(img_res->ni()-1.0, std::ceil(proj_bbox.max_x())));
104 int max_j = int(std::min(img_res->nj()-1.0, std::ceil(proj_bbox.max_y())));
105 int crop_ni = max_i - min_i + 1;
106 int crop_nj = max_j - min_j + 1;
107 if (crop_ni < 0 || crop_nj < 0) {
108 std::cerr << "Error: boxm2_roi_init_geotiff_process: the image does not contain scene!\n"
109 << "scene bbox:\n " << scene_bbox << '\n';
110 double n,e; int utm_zone;
111 std::cerr.setf(std::ios::fixed, std::ios::floatfield); std::cerr.precision(2);
112 geocam->local_to_utm(scene_bbox.min_x(), scene_bbox.min_y(), scene_bbox.min_z(), e, n, utm_zone);
113 std::cerr << "in UTM bbox transforms to " << n << " N " << e << " E and ";
114 geocam->local_to_utm(scene_bbox.max_x(), scene_bbox.max_y(), scene_bbox.max_z(), e, n, utm_zone);
115 std::cerr << n << " N " << e << " E\n";
116 return false;
117 }
118
119 vil_image_view_base_sptr crop_base = img_res->get_view((unsigned int)min_i, crop_ni, (unsigned int)min_j, crop_nj);
120 if (!crop_base) {
121 std::cerr << "Error: boxm2_roi_init_geotiff_process: could not crop the scene in this image!\n";
122 return false;
123 }
124
125 // downsample the image based on the level
126 if (level < 32) {
127 crop_ni /= (1L<<level);
128 crop_nj /= (1L<<level);
129 }
130 else {
131 crop_ni = (unsigned int) (crop_ni*std::pow(0.5,static_cast<double>(level)));
132 crop_nj = (unsigned int) (crop_nj*std::pow(0.5,static_cast<double>(level)));
133 }
134
135 // now use this camera to generate generic cam
136 vpgl_generic_camera<double> gcam;
137 // make the geocam local
138 geocam->translate(min_i, min_j, 0);
139
140 vpgl_generic_camera_convert::convert(*geocam, crop_ni, crop_nj, scene_bbox.max_z()+1.0, gcam, level);
141
142 vpgl_camera_double_sptr out = new vpgl_generic_camera<double>(gcam);
143 pro.set_output_val<vpgl_camera_double_sptr>(0, out);
144
145 auto* crop_view = dynamic_cast<vil_image_view<vxl_byte>*>(crop_base.ptr());
146 if (!crop_view) {
147 auto* crop_view_f = dynamic_cast<vil_image_view<float>*>(crop_base.ptr());
148 if (!crop_view_f) {
149 std::cerr << "Error: boxm2_roi_init_geotiff_process: could not cast first return image to a vil_image_view<float>\n";
150 return false;
151 }
152 auto* out_img = new vil_image_view<float>(crop_ni, crop_nj, crop_view_f->nplanes());
153 vil_resample_bilin(*crop_view_f, *out_img, crop_ni, crop_nj);
154 pro.set_output_val<vil_image_view_base_sptr>(1,new vil_image_view<float>(*out_img));
155 return true;
156 }
157
158 //: downsample the image based on the level
159 auto* out_img = new vil_image_view<vxl_byte>(crop_ni, crop_nj, crop_view->nplanes());
160 vil_resample_bilin(*crop_view, *out_img, crop_ni, crop_nj);
161 pro.set_output_val<vil_image_view_base_sptr>(1,new vil_image_view<vxl_byte>(*out_img));
162
163 return true;
164 }
165