Skip to content

Commit 2fcf01f

Browse files
author
Erika
committed
added calculateWeight function
1 parent e690254 commit 2fcf01f

File tree

3 files changed

+164
-60
lines changed

3 files changed

+164
-60
lines changed

‎depthmap/depthmap.cpp

+87-39
Original file line numberDiff line numberDiff line change
@@ -163,7 +163,7 @@ Eigen::Vector3f Camera::projectionToImage(Eigen::Vector3f realPosition) const{
163163

164164
//proiezione
165165
if (cameraCoords.z() == 0) {
166-
cerr << "Warning: Z è zero, impossibile proiettare il punto." << endl;
166+
cerr << "Z è zero, impossibile proiettare il punto." << endl;
167167
return Eigen::Vector3f(0, 0, 0);
168168
}
169169
//Normalize by dividing by the z coordinate to get the image coordinates u and v as projected 2D coordinates
@@ -1001,7 +1001,7 @@ void OrthoDepthmap::integratedCamera(const CameraDepthmap& camera, const char *o
10011001

10021002
for (size_t i = 0; i < point_cloud.size(); i++) {
10031003

1004-
Eigen::Vector3f& realCoord = point_cloud[i];
1004+
Eigen::Vector3f realCoord = point_cloud[i];
10051005
float h = realCoord[2];
10061006
Eigen::Vector3f pixelCoord = realToPixelCoord(realCoord[0], realCoord[1], realCoord[2]);
10071007
// project from ortho plane to camera plane, hence the fixed z
@@ -1013,17 +1013,7 @@ void OrthoDepthmap::integratedCamera(const CameraDepthmap& camera, const char *o
10131013

10141014
if (pixelX >= 0 && pixelX < camera.width && pixelY >= 0 && pixelY < camera.height) {
10151015
float depthValue = camera.elevation[pixelX + pixelY * camera.width];
1016-
if(fabs(pixelCoord[0] - 170) < 3 && fabs(pixelCoord[1] - 150) < 3){
1017-
cout << "elevation: " << elevation[int(pixelCoord[0]) + int(pixelCoord[1]) *width] << endl;
1018-
1019-
cout << "PixelX: " << pixelCoord[0] << ", PixelY: " << pixelCoord[1] << "pixelZ " << pixelCoord[2] << endl;
1020-
cout << "h: " << h << endl;
1021-
cout << "depthValue: " << depthValue << endl;
1022-
cout << "Image Coords X: " << imageCoords[0] << ", Y: " << imageCoords[1] << endl;
1023-
cout << "Camera Width: " << camera.width << ", Camera Height: " << camera.height << endl;
1024-
}
1025-
1026-
1016+
cout << pixelX << " " << pixelY << " " << depthValue << endl;
10271017
imageCloud.push_back(Eigen::Vector3f(imageCoords[0]/camera.width, imageCoords[1]/camera.height, h));
10281018
source.push_back(depthValue);
10291019
out << depthValue << "\t" << h << "\t" << imageCoords[0]/camera.width << "\t" << imageCoords[1]/camera.height <<"\n"; //red e green
@@ -1050,31 +1040,24 @@ void OrthoDepthmap::integratedCamera(const CameraDepthmap& camera, const char *o
10501040

10511041
if(px >= 0 && px< camera.width && py >= 0 && py< camera.height){
10521042
float depthValue = camera.elevation[px + py * camera.width];
1053-
p[0] /= camera.width;
1054-
p[1] /= camera.height;
1043+
p[0] = px / float(camera.width);
1044+
p[1] = py / float(camera.height);
10551045
float h = gaussianGrid.target(p[0], p[1], depthValue);
10561046
Eigen::Vector3f d = realToPixelCoord(r[0], r[1], h);
10571047

1048+
float w = camera.calculateWeight(px, py);
1049+
10581050
//p0 e p1 devono venire uguale e vedi se depth è ugusle, h dovrebbe venire simile
1059-
if (fabs(x - 170) < 1 && fabs(y - 150) < 1){
1060-
cout << "Real Coord r:" << r[0] << endl;
1061-
cout << "Projected Coords p:" << p[0] << endl;
1062-
cout << "comparison:" << endl;
1063-
cout << "PixelX: " << px << ", PixelY: " << py << endl;
1064-
cout << "p[0]: " << p[0] << ", p[1]: " << p[1] << endl;
1065-
cout << "Depth Value: " << depthValue << endl;
1066-
cout <<"gaussian grid" << gaussianGrid.value(p[0], p[1]) << endl;
1067-
cout << "elevation: " << elevation[x+y * width] << endl;
1068-
cout << " h: " << h << endl;
1069-
cout << " d: " << d << endl;
1070-
}
1071-
elevation[x + y * width] = d[2];
1051+
elevation[x + y * width] += w * d[2];
1052+
mask[x+ y * width] += w;
10721053
} else {
10731054
//elevation[x + y*width] = origin[2] + resolution[2] * elevation[x + y*width];
10741055
}
10751056
}
10761057
}
1058+
saveObj("testElev.obj");
10771059
}
1060+
/*_-----------------------------------------------------------------------------------------*/
10781061
void GaussianGrid::fitLinear(std::vector<float> &x, std::vector<float> &y, float &a, float &b) {
10791062
if (x.size() != y.size()) {
10801063
cout << "Errore: i vettori x e y devono avere la stessa lunghezza." << endl;
@@ -1096,6 +1079,32 @@ void GaussianGrid::fitLinear(std::vector<float> &x, std::vector<float> &y, float
10961079
b = (sum_y - a * sum_x) / n;
10971080
}
10981081

1082+
float GaussianGrid::bilinearInterpolation(float x, float y) {
1083+
1084+
float x1 = floor(x);
1085+
float y1 = floor(y);
1086+
float x2 = x1+1;
1087+
float y2 = y1+1;
1088+
1089+
1090+
if (x1 < 0 || x2 >= width || y1 < 0 || y2 >= height) {
1091+
cerr << "Coordinate fuori dai limiti della griglia!" << endl;
1092+
return 0.0f;
1093+
}
1094+
1095+
1096+
float Q11 = values[x1 + y1 * width];
1097+
float Q12 = values[x1 + y2 * width];
1098+
float Q21 = values[x2 + y1 * width];
1099+
float Q22 = values[x2 + y2 * width];
1100+
1101+
float R1 = (x2 - x) * Q11 + (x - x1) * Q21;
1102+
float R2 = (x2 - x) * Q12 + (x - x1) * Q22;
1103+
1104+
float P = (y2 - y) * R1 + (y - y1) * R2;
1105+
1106+
return P;
1107+
}
10991108
//fit h = a+b*elev
11001109
void GaussianGrid::init(std::vector<Eigen::Vector3f> &cloud, std::vector<float> &source) {
11011110
int side = static_cast<int>(sqrt(cloud.size())) / 2;
@@ -1113,6 +1122,7 @@ void GaussianGrid::init(std::vector<Eigen::Vector3f> &cloud, std::vector<float>
11131122
for(size_t i = 0; i < cloud.size(); i++) {
11141123
cloud[i][2] -= depthmapToCloud(source[i]);
11151124
}
1125+
11161126
computeGaussianWeightedGrid(cloud);
11171127
fillLaplacian(precision);
11181128
}
@@ -1207,14 +1217,18 @@ void GaussianGrid::fillLaplacian(float precision){
12071217
float GaussianGrid::value(float x, float y){
12081218
//bicubic interpolation
12091219
//nearest
1210-
int pixelX = std::max(0, std::min(width-1, (int)(x * (width-1))));
1211-
int pixelY = std::max(0, std::min(height-1, (int)(y * (height-1))));
1212-
return values[pixelX + pixelY * width];
1220+
float pixelX = x * (width-1);
1221+
float pixelY = y * (height-1);
1222+
1223+
return bilinearInterpolation(pixelX, pixelY);
1224+
1225+
//return values[pixelX + pixelY * width];
1226+
12131227
}
12141228

12151229
float GaussianGrid::target(float x, float y, float h) {
12161230
h = depthmapToCloud(h);
1217-
return h - value(x, y);
1231+
return h + value(x, y);
12181232
}
12191233

12201234
void GaussianGrid::computeGaussianWeightedGrid(std::vector<Eigen::Vector3f> &differences) {
@@ -1230,8 +1244,10 @@ void GaussianGrid::computeGaussianWeightedGrid(std::vector<Eigen::Vector3f> &dif
12301244
values.resize(width * height, 0);
12311245
weights.resize(width * height, 0);
12321246

1247+
std::vector<int> count(width*height, 0);
12331248

1234-
float max_distance = 3 * sigma;
1249+
1250+
float max_distance = 1.2 * sigma;
12351251
for (auto &p : differences) {
12361252

12371253

@@ -1240,22 +1256,30 @@ void GaussianGrid::computeGaussianWeightedGrid(std::vector<Eigen::Vector3f> &dif
12401256
int y_start = std::max(0, static_cast<int>((p[1] - max_distance - y_min) / y_step));
12411257
int y_end = std::min(height - 1, static_cast<int>((p[1] + max_distance - y_min) / y_step));
12421258

1243-
for (int i = x_start; i <= x_end; i++) {
1244-
for (int j = y_start; j <= y_end; j++) {
1245-
float xg = x_min + i * x_step;
1246-
float yg = y_min + j * y_step;
1259+
for (int x = x_start; x <= x_end; x++) {
1260+
for (int y = y_start; y <= y_end; y++) {
1261+
float xg = x_min + x * x_step;
1262+
float yg = y_min + y * y_step;
12471263

12481264
float distance = sqrt((p[0] - xg) * (p[0] - xg) + (p[1] - yg) * (p[1] - yg));
12491265
if (distance <= max_distance) {
12501266
float weight = exp(-(distance * distance) / (2 * sigma * sigma));
1251-
values[i * width + j] += weight * p[2];
1252-
weights[i * width + j] += weight;
1267+
values[y * width + x] += weight * p[2];
1268+
weights[y * width + x] += weight;
1269+
count[y*width + x]++;
12531270
}
12541271
}
12551272
}
12561273
}
1274+
//chiama camere tutte e 4 e vedi come vengono
1275+
// pesare per il blanding funzione intervallo 0, 1 * 0, 1 0 ai bordi 1 al centro, che sia una funzione continua
1276+
//polinomio di 2 grado in x * pol 2 grado in y. derivata e peso a 0 sul bordo
1277+
//fai somma pesata e veedi come vieni
1278+
//funz target ritorna valore e peso
12571279

12581280
for (int i = 0; i < values.size(); i++) {
1281+
if(count[i] < 3)
1282+
weights[i] = 0;
12591283
if (weights[i] != 0) {
12601284
values[i] /= (weights[i]);
12611285
}
@@ -1266,7 +1290,31 @@ void GaussianGrid::computeGaussianWeightedGrid(std::vector<Eigen::Vector3f> &dif
12661290
}
12671291
}*/
12681292
}
1293+
//
1294+
float Depthmap::calculateWeight(float x, float y) const{
1295+
1296+
x/=width;
1297+
y/= height;
1298+
1299+
float weightX =pow(cos(M_PI * (x-0.5f)), 2);
1300+
float weightY = pow(cos(M_PI * (y-0.5f)), 2);
1301+
1302+
return weightX * weightY;
1303+
}
1304+
void OrthoDepthmap::beginIntegration(){
1305+
elevation.clear();
1306+
elevation.resize(width * height, 0);
1307+
mask.clear();
1308+
mask.resize(width * height, 0);
12691309

1310+
}
1311+
void OrthoDepthmap::endIntegration(){
1312+
for(size_t i =0; i < elevation.size(); i++){
1313+
if(mask[i]){
1314+
elevation[i] /= mask[i];
1315+
}
1316+
}
1317+
}
12701318

12711319
void GaussianGrid::imageGrid(const char* filename) {
12721320

‎depthmap/depthmap.h

+6-1
Original file line numberDiff line numberDiff line change
@@ -50,6 +50,7 @@ class Depthmap {
5050
void computeNormals();
5151
void depthIntegrateNormals();
5252
void resizeNormals(int factorPowerOfTwo, int step = 1);
53+
float calculateWeight(float x, float y) const; //
5354

5455

5556
protected:
@@ -70,14 +71,16 @@ class GaussianGrid {
7071

7172
void init(std::vector<Eigen::Vector3f> &cloud, std::vector<float> &source);
7273
void fitLinear(std::vector<float> &x, std::vector<float> &y, float &a, float &b); //ax + b
73-
74+
float bilinearInterpolation(float x, float y);
7475
void computeGaussianWeightedGrid(std::vector<Eigen::Vector3f> &differences);
7576
void fillLaplacian(float precision);
7677
void imageGrid(const char* filename);
7778
//interpola la griglia per spostare la depthmap, serve per creare la funzione
7879
float value(float x, float y);
7980
float target(float x, float y, float source); //given a point in the source depthmap compute the z in cloud coordinate space;
8081

82+
83+
8184
float depthmapToCloud(float h) {
8285
return a*h + b;
8386
}
@@ -117,6 +120,8 @@ class OrthoDepthmap:
117120
//legge nella depth l h corrispondente
118121
void verifyPointCloud();
119122
void integratedCamera(const CameraDepthmap& camera, const char *outputFile);
123+
void beginIntegration();
124+
void endIntegration();
120125

121126

122127
/*1.

‎depthmap/main.cpp

+71-20
Original file line numberDiff line numberDiff line change
@@ -2,6 +2,8 @@
22
#include "../src/bni_normal_integration.h"
33
#include <iostream>
44
#include <QFile>
5+
#include <QDir>
6+
#include <QFileInfoList>
57
#include <QDomDocument>
68
#include <eigen3/Eigen/Dense>
79
#include <math.h>
@@ -15,19 +17,21 @@ int main(int argc, char *argv[]) {
1517
return 1;
1618
}*/
1719
//input
18-
//#define MACOS 1
20+
#define MACOS 1
1921
#ifdef MACOS
20-
QString base = "/Users/erika/Desktop/";
22+
QString base = "/Users/erika/Desktop/testcenterRel_copia/";
2123
#else
2224
QString base = "";
2325
#endif
2426

25-
QString depthmapPath = base + "Z_Num7_DeZoom4_STD-MALT.tif";
26-
QString cameraDepthmap = base + "L04C12_depth_rti.tiff";
27-
QString orientationXmlPath = base + "Orientation-L04C12.tif.xml";
28-
QString maskPath = base + "Masq_STD-MALT_DeZoom4.tif";
29-
QString plyFile = base +"AperiCloud_Relative__mini.ply";
30-
QString point_txt = base + "points_h.txt";
27+
28+
29+
QString depthmapPath = base + "photogrammetry/Malt/Z_Num7_DeZoom4_STD-MALT.tif";
30+
//QString cameraDepthmap = base + "datasets/L04C12.tif";
31+
//QString orientationXmlPath = base + "photogrammetry/Ori-Relative/Orientation-L04C12.tif.xml";
32+
QString maskPath = base + "photogrammetry/Malt/Masq_STD-MALT_DeZoom4.tif";
33+
QString plyFile = base +"photogrammetry/AperiCloud_Relative__mini.ply";
34+
QString point_txt = base + "photogrammetry/points_h.txt";
3135
Depthmap depth;
3236

3337
//output
@@ -43,6 +47,24 @@ int main(int argc, char *argv[]) {
4347
cout << "accidenti" << endl;
4448
return -1;
4549
}
50+
51+
QDir datasetsDir(base + "datasets");
52+
QDir xmlDir(base + "photogrammetry/Ori-Relative");
53+
QStringList tiffFilters = {"*.tif"};
54+
QStringList xmlFilters = {"Orientation-*.tif.xml"};
55+
56+
QFileInfoList tiffFiles = datasetsDir.entryInfoList(tiffFilters, QDir::Files);
57+
if (tiffFiles.isEmpty()) {
58+
cerr << "No .tiff files found in " << datasetsDir.absolutePath().toStdString() << endl;
59+
return -1;
60+
}
61+
62+
QFileInfoList xmlFiles = xmlDir.entryInfoList(xmlFilters, QDir::Files);
63+
if (xmlFiles.isEmpty()) {
64+
cerr << "No .xml files found in " << xmlDir.absolutePath().toStdString() << endl;
65+
return -1;
66+
}
67+
4668
//ortho.computeNormals();
4769
//ortho.saveNormals(qPrintable(base + "testcenterRel_copia/photogrammetry/original.png"));
4870
//ortho.saveObj(qPrintable(base + "testcenterRel_copia/photogrammetry/original.obj"));
@@ -54,37 +76,66 @@ int main(int argc, char *argv[]) {
5476
cout << qPrintable(e) << endl;
5577
exit(-1);
5678
}
79+
5780
ortho.verifyPointCloud();
81+
ortho.beginIntegration();
82+
83+
84+
for (const QFileInfo &tiffFile : tiffFiles) {
85+
86+
CameraDepthmap depthCam;
87+
QString cameraName = tiffFile.completeBaseName();
88+
QString orientationXmlPath = xmlDir.absoluteFilePath("Orientation-" + cameraName + ".tif.xml");
89+
90+
cout << "Looking for XML: " << orientationXmlPath.toStdString() << endl;
91+
92+
if (!depthCam.camera.loadXml(orientationXmlPath)) {
93+
cerr << "Failed to load XML: " << orientationXmlPath.toStdString() << endl;
94+
continue;
95+
}
96+
97+
if (!depthCam.loadDepth(qPrintable(tiffFile.absoluteFilePath()))) {
98+
cerr << "Failed to load depth map: " << tiffFile.fileName().toStdString() << endl;
99+
continue;
100+
}
101+
if(depthCam.width != depthCam.camera.width || depthCam.height != depthCam.camera.height){
102+
cerr << "width is not the same" << endl;
103+
return -1;
104+
}
105+
cout << "Processed: " << tiffFile.fileName().toStdString() << endl;
106+
107+
ortho.integratedCamera(depthCam, qPrintable(output_points));
108+
//ortho.projectToCameraDepthMap(depthCam.camera, outputPath);
109+
110+
QString outputTiffPath = base +"output_" + tiffFile.fileName();
111+
cout << "Output TIFF Path: " << outputTiffPath.toStdString() << endl;
112+
58113

59-
CameraDepthmap depthCam;
60114

61-
//xml per camera e tiff per la depth map
62-
depthCam.camera.loadXml(orientationXmlPath);
63-
depthCam.loadDepth(qPrintable(cameraDepthmap));
64-
if(depthCam.width != depthCam.camera.width || depthCam.height != depthCam.camera.height){
65-
cerr << "width is not the same" << endl;
66-
return -1;
67115
}
68-
ortho.integratedCamera(depthCam, qPrintable(output_points));
69-
ortho.saveDepth(qPrintable(base + "testDepth.tiff"));
70-
// sqrt
116+
ortho.endIntegration();
117+
ortho.saveDepth(qPrintable("final_depth.tif"));
71118

72119

120+
//depthCam.camera.loadXml(orientationXmlPath);
121+
//depthCam.loadDepth(qPrintable(tiffFile.absoluteFilePath()));
122+
123+
//ortho.saveDepth(qPrintable(base + "testDepth.tiff"));
73124
//ortho.computeGaussianWeightedGrid(qPrintable(point_txt));
74125

75126

76127
//int pixelX = 165;
77128
//int pixelY = 144;
78129
//float pixelZ = 4.5;
79-
ortho.projectToCameraDepthMap(depthCam.camera, outputPath);
130+
80131
Eigen::Matrix3f rotationMatrix;
81132
Eigen::Vector3f center;
82133

83134
// depth.loadMask(qPrintable(maskPath));
84135
// depth.saveDepth(qPrintable(depthmapPath));
85136
// depth.saveMask(qPrintable(maskPath));
86137
//QString maskObjPath = base + "testcenterRel_copia/photogrammetry/mask.obj";
87-
ortho.saveObj(qPrintable(base + "testcenterRel_copia/photogrammetry/depthmap_projectL05C13.obj"));
138+
ortho.saveObj(qPrintable(base + "depthmap_projectL05C13.obj"));
88139

89140

90141
//depth.depthIntegrateNormals();

0 commit comments

Comments
 (0)