From 7a436949b56d1bb166feaa4c2944645051f8049b Mon Sep 17 00:00:00 2001 From: ortrr Date: Thu, 23 Apr 2026 19:04:35 +0300 Subject: [PATCH] Done: task05 --- src/phg/mvs/depth_maps/pm_depth_maps.cpp | 90 +++++++++++++++--------- tests/test_depth_maps_pm.cpp | 68 +++++++++--------- 2 files changed, 92 insertions(+), 66 deletions(-) diff --git a/src/phg/mvs/depth_maps/pm_depth_maps.cpp b/src/phg/mvs/depth_maps/pm_depth_maps.cpp index 7021ddd0..0bbd1b08 100644 --- a/src/phg/mvs/depth_maps/pm_depth_maps.cpp +++ b/src/phg/mvs/depth_maps/pm_depth_maps.cpp @@ -46,18 +46,20 @@ vector3d project(const vector3d& global_point, const phg::Calibration& calibrati return pixel_with_depth; } -// TODO 101 реализуйте unproject (вам поможет тест на идемпотентность project -> unproject в test_depth_maps_pm) -vector3d unproject(const vector3d& pixel, const phg::Calibration& calibration, const matrix34d& PtoWorld) +vector3d unproject(const vector3d& pix, const phg::Calibration& calib, const matrix34d& worldFromPose) { - double depth = pixel[2]; // на самом деле это не глубина, это координата по оси +Z (вдоль которой смотрит камера в ее локальной системе координат) + double z_cam = pix[2]; // локальная Z-координата (вдоль оптической оси камеры) - vector3d local_point; // TODO 102 пустите луч pixel из calibration а затем возьмите ан нем точку у которой по оси +Z координата=depth + vector3d tmp = calib.unproject({ pix[0], pix[1] }); - vector3d global_point; // TODO 103 переведите точку из локальной системы в глобальную + const double factor = z_cam / tmp[2]; + vector3d local = tmp * factor; - return global_point; + vector3d global = worldFromPose * homogenize(local); + return global; } + void PMDepthMapsBuilder::buildDepthMap(unsigned int camera_key, cv::Mat& depth_map_res, cv::Mat& normal_map_res, cv::Mat& cost_map_res, float depth_min, float depth_max) { rassert(camera_key < ncameras, 238192841294108); @@ -112,16 +114,18 @@ void PMDepthMapsBuilder::refinement() n0 = normal_map.at(j, i); // 2) случайной пертурбации текущей гипотезы (мутация и уточнение того что уже смогли найти) - dp = r.nextf(d0 * 0.5f, d0 * 1.5); // TODO 104: сделайте так чтобы отклонение было тем меньше, чем номер итерации ближе к NITERATIONS, улучшило ли это результат? - np = cv::normalize(n0 + randomNormalObservedFromCamera(cameras_RtoWorld[ref_cam], r) * 0.5); // TODO 105: сделайте так чтобы отклонение было тем меньше, чем номер итерации ближе к NITERATIONS, улучшило ли это результат? + float t = float(iter) / NITERATIONS; // t = 0..1, где 1 — финальная итерация + float min_factor = std::max(0.1f, 0.5f * (1.0f - t)); + float max_factor = std::max(0.2f, 1.5f * (1.0f - t)); + dp = r.nextf(d0 * min_factor, d0 * max_factor); + + np = cv::normalize(n0 + randomNormalObservedFromCamera(cameras_RtoWorld[ref_cam], r) * 0.5); dp = std::max(ref_depth_min, std::min(ref_depth_max, dp)); // 3) новой случайной гипотезы из фрустума поиска (новые идеи, вечный поиск во всем пространстве) - // TODO 106: создайте случайную гипотезу dr+nr, вам поможет: - // - r.nextf(...) - // - ref_depth_min, ref_depth_max - // - randomNormalObservedFromCamera - поможет создать нормаль которая гарантированно смотрит на нас + dr = r.nextf(ref_depth_min, ref_depth_max); + nr = randomNormalObservedFromCamera(cameras_RtoWorld[ref_cam], r); } float best_depth = d0; @@ -315,10 +319,13 @@ float PMDepthMapsBuilder::estimateCost(ptrdiff_t i, ptrdiff_t j, double d, const vector3d pixel(i + 0.5, j + 0.5, d); vector3d global_point = unproject(pixel, calibration, cameras_PtoWorld[ref_cam]); - if (!(i - COST_PATCH_RADIUS >= 0 && i + COST_PATCH_RADIUS < width)) + if (!(i - COST_PATCH_RADIUS >= 0 && i + COST_PATCH_RADIUS < width)) { return NO_COST; - if (!(j - COST_PATCH_RADIUS >= 0 && j + COST_PATCH_RADIUS < height)) + } + + if (!(j - COST_PATCH_RADIUS >= 0 && j + COST_PATCH_RADIUS < height)) { return NO_COST; + } std::vector patch0, patch1; @@ -348,35 +355,48 @@ float PMDepthMapsBuilder::estimateCost(ptrdiff_t i, ptrdiff_t j, double d, const double x = neighb_proj[0]; double y = neighb_proj[1]; - // TODO 205: замените этот наивный вариант nearest neighbor сэмплирования текстуры на билинейную интерполяцию (учтите что центр пикселя - .5 после запятой) - ptrdiff_t u = x; - ptrdiff_t v = y; + ptrdiff_t u = static_cast(x); + ptrdiff_t v = static_cast(y); - // TODO 108: добавьте проверку "попали ли мы в камеру номер neighb_cam?" если не попали - возвращаем NO_COST + if (u < 0 || u >= width || v < 0 || v >= height) { + return NO_COST; + } float intensity = cameras_imgs_grey[neighb_cam].at(v, u) / 255.0f; patch1.push_back(intensity); } } - // TODO 109: реализуйте ZNCC https://en.wikipedia.org/wiki/Cross-correlation#Zero-normalized_cross-correlation_(ZNCC) - // или слайд #25 в лекции 5 про SGM и Cost-функции - https://my.compscicenter.ru/attachments/classes/slides_w2n8WNLY/photogrammetry_lecture_090321.pdf + // ZNCC rassert(patch0.size() == patch1.size(), 12489185129326); size_t n = patch0.size(); float mean0 = 0.0f; float mean1 = 0.0f; - // ... + for (size_t k = 0; k < n; ++k) { - float a = patch0[k]; - float b = patch1[k]; - mean0 += a; - mean1 += b; - // ... + mean0 += patch0[k]; + mean1 += patch1[k]; } mean0 /= n; mean1 /= n; - // ... + + float sum_prod = 0.0f; + float sum_sq0 = 0.0f; + float sum_sq1 = 0.0f; + + for (size_t k = 0; k < n; ++k) { + float a = patch0[k] - mean0; + float b = patch1[k] - mean1; + sum_prod += a * b; + sum_sq0 += a * a; + sum_sq1 += b * b; + } + float zncc = 0.0f; + float denom = std::sqrt(sum_sq0 * sum_sq1); + if (denom > 1e-6f) { + zncc = sum_prod / denom; + } // ZNCC в диапазоне [-1; 1], 1: идеальное совпадение, -1: ничего общего rassert(zncc == zncc, 23141241210380); // проверяем что не nan @@ -393,8 +413,9 @@ float PMDepthMapsBuilder::estimateCost(ptrdiff_t i, ptrdiff_t j, double d, const float PMDepthMapsBuilder::avgCost(std::vector& costs) { - if (costs.size() == 0) + if (costs.size() == 0) { return NO_COST; + } std::sort(costs.begin(), costs.end()); @@ -403,10 +424,15 @@ float PMDepthMapsBuilder::avgCost(std::vector& costs) float cost_sum = best_cost; float cost_w = 1.0f; - // TODO 110 реализуйте какое-то "усреднение cost-ов по всем соседям", с ограничением что участвуют только COSTS_BEST_K_LIMIT лучших - // TODO 111 добавьте к этому усреднению еще одно ограничение: если cost больше чем best_cost*COSTS_K_RATIO - то такой cost подозрительно плохой и мы его не хотим учитывать (вероятно occlusion) - // TODO 112 а что если в пикселе occlusion, но best_cost - большой и поэтому отсечение по best_cost*COSTS_K_RATIO не срабатывает? можно ли это отсечение как-то выправить для такого случая? - // TODO 207 а что если добавить какой-нибудь бонус в случае если больше чем Х камер засчиталось? улучшается/ухудшается ли от этого что-то на herzjezu25? а при большем числе фотографий + const int k_limit = std::min((int)costs.size(), COSTS_BEST_K_LIMIT); + for (int i = 1; i < k_limit; ++i) { + float c = costs[i]; + if (c > best_cost * COSTS_K_RATIO) { + break; // последующие cost'ы только больше (так как массив отсортирован) + } + cost_sum += c; + cost_w += 1.0f; + } float avg_cost = cost_sum / cost_w; return avg_cost; diff --git a/tests/test_depth_maps_pm.cpp b/tests/test_depth_maps_pm.cpp index 3a03547e..e64226d6 100644 --- a/tests/test_depth_maps_pm.cpp +++ b/tests/test_depth_maps_pm.cpp @@ -34,43 +34,43 @@ TEST(test_depth_maps_pm, SingleDepthMap) { // TODO этот код надо раскомментировать чтобы запустить тестирование: - // Dataset dataset = loadDataset(DATASET_DIR, DATASET_DOWNSCALE); - // phg::PMDepthMapsBuilder builder(dataset.ncameras, dataset.cameras_imgs, dataset.cameras_imgs_grey, dataset.cameras_labels, dataset.cameras_P, dataset.calibration); - // - // size_t ci = 2; // строим карту глубины для третьей камеры (индексация с нуля) - // size_t cameras_limit = 5; // учитывая первые пять фотографий датасета, т.е. две камеры слева и две камеры справа - // - // dataset.ncameras = cameras_limit; - // cv::Mat depth_map, normal_map, cost_map; - // builder.buildDepthMap(ci, depth_map, cost_map, normal_map, dataset.cameras_depth_min[ci], dataset.cameras_depth_max[ci]); + Dataset dataset = loadDataset(DATASET_DIR, DATASET_DOWNSCALE); + phg::PMDepthMapsBuilder builder(dataset.ncameras, dataset.cameras_imgs, dataset.cameras_imgs_grey, dataset.cameras_labels, dataset.cameras_P, dataset.calibration); + + size_t ci = 2; // строим карту глубины для третьей камеры (индексация с нуля) + size_t cameras_limit = 5; // учитывая первые пять фотографий датасета, т.е. две камеры слева и две камеры справа + + dataset.ncameras = cameras_limit; + cv::Mat depth_map, normal_map, cost_map; + builder.buildDepthMap(ci, depth_map, cost_map, normal_map, dataset.cameras_depth_min[ci], dataset.cameras_depth_max[ci]); } TEST(test_depth_maps_pm, AllDepthMaps) { // TODO этот код можно раскомментировать чтобы построить много карт глубины и сохранить их облака точек: - // Dataset full_dataset = loadDataset(DATASET_DIR, DATASET_DOWNSCALE); - // - // const size_t ref_camera_shift = 2; - // const size_t to_shift = 5; - // - // std::vector all_points; - // std::vector all_colors; - // std::vector all_normals; - // - // size_t ndepth_maps = 0; - // - // for (size_t from = 0; from + to_shift <= full_dataset.ncameras; ++from) { - // size_t to = from + to_shift; - // - // Dataset dataset = full_dataset.subset(from, to); - // - // phg::PMDepthMapsBuilder builder(dataset.ncameras, dataset.cameras_imgs, dataset.cameras_imgs_grey, dataset.cameras_labels, dataset.cameras_P, dataset.calibration); - // cv::Mat depth_map, normal_map, cost_map; - // builder.buildDepthMap(ref_camera_shift, depth_map, cost_map, normal_map, dataset.cameras_depth_min[ref_camera_shift], dataset.cameras_depth_max[ref_camera_shift]); - // phg::PMDepthMapsBuilder::buildGoodPoints(depth_map, normal_map, cost_map, dataset.cameras_imgs[ref_camera_shift], dataset.calibration, builder.getCameraPtoWorld(ref_camera_shift), all_points, all_colors, all_normals); - // ++ndepth_maps; - // - // std::string tie_points_filename = std::string("data/debug/test_depth_maps_pm/") + getTestName() + "/all_points_" + to_string(ndepth_maps) + ".ply"; - // phg::exportPointCloud(all_points, tie_points_filename, all_colors, all_normals); - // } + Dataset full_dataset = loadDataset(DATASET_DIR, DATASET_DOWNSCALE); + + const size_t ref_camera_shift = 2; + const size_t to_shift = 5; + + std::vector all_points; + std::vector all_colors; + std::vector all_normals; + + size_t ndepth_maps = 0; + + for (size_t from = 0; from + to_shift <= full_dataset.ncameras; ++from) { + size_t to = from + to_shift; + + Dataset dataset = full_dataset.subset(from, to); + + phg::PMDepthMapsBuilder builder(dataset.ncameras, dataset.cameras_imgs, dataset.cameras_imgs_grey, dataset.cameras_labels, dataset.cameras_P, dataset.calibration); + cv::Mat depth_map, normal_map, cost_map; + builder.buildDepthMap(ref_camera_shift, depth_map, cost_map, normal_map, dataset.cameras_depth_min[ref_camera_shift], dataset.cameras_depth_max[ref_camera_shift]); + phg::PMDepthMapsBuilder::buildGoodPoints(depth_map, normal_map, cost_map, dataset.cameras_imgs[ref_camera_shift], dataset.calibration, builder.getCameraPtoWorld(ref_camera_shift), all_points, all_colors, all_normals); + ++ndepth_maps; + + std::string tie_points_filename = std::string("data/debug/test_depth_maps_pm/") + getTestName() + "/all_points_" + to_string(ndepth_maps) + ".ply"; + phg::exportPointCloud(all_points, tie_points_filename, all_colors, all_normals); + } }