Skip to content

Commit

Permalink
added calculateWeight function
Browse files Browse the repository at this point in the history
  • Loading branch information
Erika committed Jan 22, 2025
1 parent e690254 commit 2fcf01f
Show file tree
Hide file tree
Showing 3 changed files with 164 additions and 60 deletions.
126 changes: 87 additions & 39 deletions depthmap/depthmap.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -163,7 +163,7 @@ Eigen::Vector3f Camera::projectionToImage(Eigen::Vector3f realPosition) const{

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

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

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

if (pixelX >= 0 && pixelX < camera.width && pixelY >= 0 && pixelY < camera.height) {
float depthValue = camera.elevation[pixelX + pixelY * camera.width];
if(fabs(pixelCoord[0] - 170) < 3 && fabs(pixelCoord[1] - 150) < 3){
cout << "elevation: " << elevation[int(pixelCoord[0]) + int(pixelCoord[1]) *width] << endl;

cout << "PixelX: " << pixelCoord[0] << ", PixelY: " << pixelCoord[1] << "pixelZ " << pixelCoord[2] << endl;
cout << "h: " << h << endl;
cout << "depthValue: " << depthValue << endl;
cout << "Image Coords X: " << imageCoords[0] << ", Y: " << imageCoords[1] << endl;
cout << "Camera Width: " << camera.width << ", Camera Height: " << camera.height << endl;
}


cout << pixelX << " " << pixelY << " " << depthValue << endl;
imageCloud.push_back(Eigen::Vector3f(imageCoords[0]/camera.width, imageCoords[1]/camera.height, h));
source.push_back(depthValue);
out << depthValue << "\t" << h << "\t" << imageCoords[0]/camera.width << "\t" << imageCoords[1]/camera.height <<"\n"; //red e green
Expand All @@ -1050,31 +1040,24 @@ void OrthoDepthmap::integratedCamera(const CameraDepthmap& camera, const char *o

if(px >= 0 && px< camera.width && py >= 0 && py< camera.height){
float depthValue = camera.elevation[px + py * camera.width];
p[0] /= camera.width;
p[1] /= camera.height;
p[0] = px / float(camera.width);
p[1] = py / float(camera.height);
float h = gaussianGrid.target(p[0], p[1], depthValue);
Eigen::Vector3f d = realToPixelCoord(r[0], r[1], h);

float w = camera.calculateWeight(px, py);

//p0 e p1 devono venire uguale e vedi se depth è ugusle, h dovrebbe venire simile
if (fabs(x - 170) < 1 && fabs(y - 150) < 1){
cout << "Real Coord r:" << r[0] << endl;
cout << "Projected Coords p:" << p[0] << endl;
cout << "comparison:" << endl;
cout << "PixelX: " << px << ", PixelY: " << py << endl;
cout << "p[0]: " << p[0] << ", p[1]: " << p[1] << endl;
cout << "Depth Value: " << depthValue << endl;
cout <<"gaussian grid" << gaussianGrid.value(p[0], p[1]) << endl;
cout << "elevation: " << elevation[x+y * width] << endl;
cout << " h: " << h << endl;
cout << " d: " << d << endl;
}
elevation[x + y * width] = d[2];
elevation[x + y * width] += w * d[2];
mask[x+ y * width] += w;
} else {
//elevation[x + y*width] = origin[2] + resolution[2] * elevation[x + y*width];
}
}
}
saveObj("testElev.obj");
}
/*_-----------------------------------------------------------------------------------------*/
void GaussianGrid::fitLinear(std::vector<float> &x, std::vector<float> &y, float &a, float &b) {
if (x.size() != y.size()) {
cout << "Errore: i vettori x e y devono avere la stessa lunghezza." << endl;
Expand All @@ -1096,6 +1079,32 @@ void GaussianGrid::fitLinear(std::vector<float> &x, std::vector<float> &y, float
b = (sum_y - a * sum_x) / n;
}

float GaussianGrid::bilinearInterpolation(float x, float y) {

float x1 = floor(x);
float y1 = floor(y);
float x2 = x1+1;
float y2 = y1+1;


if (x1 < 0 || x2 >= width || y1 < 0 || y2 >= height) {
cerr << "Coordinate fuori dai limiti della griglia!" << endl;
return 0.0f;
}


float Q11 = values[x1 + y1 * width];
float Q12 = values[x1 + y2 * width];
float Q21 = values[x2 + y1 * width];
float Q22 = values[x2 + y2 * width];

float R1 = (x2 - x) * Q11 + (x - x1) * Q21;
float R2 = (x2 - x) * Q12 + (x - x1) * Q22;

float P = (y2 - y) * R1 + (y - y1) * R2;

return P;
}
//fit h = a+b*elev
void GaussianGrid::init(std::vector<Eigen::Vector3f> &cloud, std::vector<float> &source) {
int side = static_cast<int>(sqrt(cloud.size())) / 2;
Expand All @@ -1113,6 +1122,7 @@ void GaussianGrid::init(std::vector<Eigen::Vector3f> &cloud, std::vector<float>
for(size_t i = 0; i < cloud.size(); i++) {
cloud[i][2] -= depthmapToCloud(source[i]);
}

computeGaussianWeightedGrid(cloud);
fillLaplacian(precision);
}
Expand Down Expand Up @@ -1207,14 +1217,18 @@ void GaussianGrid::fillLaplacian(float precision){
float GaussianGrid::value(float x, float y){
//bicubic interpolation
//nearest
int pixelX = std::max(0, std::min(width-1, (int)(x * (width-1))));
int pixelY = std::max(0, std::min(height-1, (int)(y * (height-1))));
return values[pixelX + pixelY * width];
float pixelX = x * (width-1);
float pixelY = y * (height-1);

return bilinearInterpolation(pixelX, pixelY);

//return values[pixelX + pixelY * width];

}

float GaussianGrid::target(float x, float y, float h) {
h = depthmapToCloud(h);
return h - value(x, y);
return h + value(x, y);
}

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

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

float max_distance = 3 * sigma;

float max_distance = 1.2 * sigma;
for (auto &p : differences) {


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

for (int i = x_start; i <= x_end; i++) {
for (int j = y_start; j <= y_end; j++) {
float xg = x_min + i * x_step;
float yg = y_min + j * y_step;
for (int x = x_start; x <= x_end; x++) {
for (int y = y_start; y <= y_end; y++) {
float xg = x_min + x * x_step;
float yg = y_min + y * y_step;

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

for (int i = 0; i < values.size(); i++) {
if(count[i] < 3)
weights[i] = 0;
if (weights[i] != 0) {
values[i] /= (weights[i]);
}
Expand All @@ -1266,7 +1290,31 @@ void GaussianGrid::computeGaussianWeightedGrid(std::vector<Eigen::Vector3f> &dif
}
}*/
}
//
float Depthmap::calculateWeight(float x, float y) const{

x/=width;
y/= height;

float weightX =pow(cos(M_PI * (x-0.5f)), 2);
float weightY = pow(cos(M_PI * (y-0.5f)), 2);

return weightX * weightY;
}
void OrthoDepthmap::beginIntegration(){
elevation.clear();
elevation.resize(width * height, 0);
mask.clear();
mask.resize(width * height, 0);

}
void OrthoDepthmap::endIntegration(){
for(size_t i =0; i < elevation.size(); i++){
if(mask[i]){
elevation[i] /= mask[i];
}
}
}

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

Expand Down
7 changes: 6 additions & 1 deletion depthmap/depthmap.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,6 +50,7 @@ class Depthmap {
void computeNormals();
void depthIntegrateNormals();
void resizeNormals(int factorPowerOfTwo, int step = 1);
float calculateWeight(float x, float y) const; //


protected:
Expand All @@ -70,14 +71,16 @@ class GaussianGrid {

void init(std::vector<Eigen::Vector3f> &cloud, std::vector<float> &source);
void fitLinear(std::vector<float> &x, std::vector<float> &y, float &a, float &b); //ax + b

float bilinearInterpolation(float x, float y);
void computeGaussianWeightedGrid(std::vector<Eigen::Vector3f> &differences);
void fillLaplacian(float precision);
void imageGrid(const char* filename);
//interpola la griglia per spostare la depthmap, serve per creare la funzione
float value(float x, float y);
float target(float x, float y, float source); //given a point in the source depthmap compute the z in cloud coordinate space;



float depthmapToCloud(float h) {
return a*h + b;
}
Expand Down Expand Up @@ -117,6 +120,8 @@ class OrthoDepthmap:
//legge nella depth l h corrispondente
void verifyPointCloud();
void integratedCamera(const CameraDepthmap& camera, const char *outputFile);
void beginIntegration();
void endIntegration();


/*1.
Expand Down
91 changes: 71 additions & 20 deletions depthmap/main.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,8 @@
#include "../src/bni_normal_integration.h"
#include <iostream>
#include <QFile>
#include <QDir>
#include <QFileInfoList>
#include <QDomDocument>
#include <eigen3/Eigen/Dense>
#include <math.h>
Expand All @@ -15,19 +17,21 @@ int main(int argc, char *argv[]) {
return 1;
}*/
//input
//#define MACOS 1
#define MACOS 1
#ifdef MACOS
QString base = "/Users/erika/Desktop/";
QString base = "/Users/erika/Desktop/testcenterRel_copia/";
#else
QString base = "";
#endif

QString depthmapPath = base + "Z_Num7_DeZoom4_STD-MALT.tif";
QString cameraDepthmap = base + "L04C12_depth_rti.tiff";
QString orientationXmlPath = base + "Orientation-L04C12.tif.xml";
QString maskPath = base + "Masq_STD-MALT_DeZoom4.tif";
QString plyFile = base +"AperiCloud_Relative__mini.ply";
QString point_txt = base + "points_h.txt";


QString depthmapPath = base + "photogrammetry/Malt/Z_Num7_DeZoom4_STD-MALT.tif";
//QString cameraDepthmap = base + "datasets/L04C12.tif";
//QString orientationXmlPath = base + "photogrammetry/Ori-Relative/Orientation-L04C12.tif.xml";
QString maskPath = base + "photogrammetry/Malt/Masq_STD-MALT_DeZoom4.tif";
QString plyFile = base +"photogrammetry/AperiCloud_Relative__mini.ply";
QString point_txt = base + "photogrammetry/points_h.txt";
Depthmap depth;

//output
Expand All @@ -43,6 +47,24 @@ int main(int argc, char *argv[]) {
cout << "accidenti" << endl;
return -1;
}

QDir datasetsDir(base + "datasets");
QDir xmlDir(base + "photogrammetry/Ori-Relative");
QStringList tiffFilters = {"*.tif"};
QStringList xmlFilters = {"Orientation-*.tif.xml"};

QFileInfoList tiffFiles = datasetsDir.entryInfoList(tiffFilters, QDir::Files);
if (tiffFiles.isEmpty()) {
cerr << "No .tiff files found in " << datasetsDir.absolutePath().toStdString() << endl;
return -1;
}

QFileInfoList xmlFiles = xmlDir.entryInfoList(xmlFilters, QDir::Files);
if (xmlFiles.isEmpty()) {
cerr << "No .xml files found in " << xmlDir.absolutePath().toStdString() << endl;
return -1;
}

//ortho.computeNormals();
//ortho.saveNormals(qPrintable(base + "testcenterRel_copia/photogrammetry/original.png"));
//ortho.saveObj(qPrintable(base + "testcenterRel_copia/photogrammetry/original.obj"));
Expand All @@ -54,37 +76,66 @@ int main(int argc, char *argv[]) {
cout << qPrintable(e) << endl;
exit(-1);
}

ortho.verifyPointCloud();
ortho.beginIntegration();


for (const QFileInfo &tiffFile : tiffFiles) {

CameraDepthmap depthCam;
QString cameraName = tiffFile.completeBaseName();
QString orientationXmlPath = xmlDir.absoluteFilePath("Orientation-" + cameraName + ".tif.xml");

cout << "Looking for XML: " << orientationXmlPath.toStdString() << endl;

if (!depthCam.camera.loadXml(orientationXmlPath)) {
cerr << "Failed to load XML: " << orientationXmlPath.toStdString() << endl;
continue;
}

if (!depthCam.loadDepth(qPrintable(tiffFile.absoluteFilePath()))) {
cerr << "Failed to load depth map: " << tiffFile.fileName().toStdString() << endl;
continue;
}
if(depthCam.width != depthCam.camera.width || depthCam.height != depthCam.camera.height){
cerr << "width is not the same" << endl;
return -1;
}
cout << "Processed: " << tiffFile.fileName().toStdString() << endl;

ortho.integratedCamera(depthCam, qPrintable(output_points));
//ortho.projectToCameraDepthMap(depthCam.camera, outputPath);

QString outputTiffPath = base +"output_" + tiffFile.fileName();
cout << "Output TIFF Path: " << outputTiffPath.toStdString() << endl;


CameraDepthmap depthCam;

//xml per camera e tiff per la depth map
depthCam.camera.loadXml(orientationXmlPath);
depthCam.loadDepth(qPrintable(cameraDepthmap));
if(depthCam.width != depthCam.camera.width || depthCam.height != depthCam.camera.height){
cerr << "width is not the same" << endl;
return -1;
}
ortho.integratedCamera(depthCam, qPrintable(output_points));
ortho.saveDepth(qPrintable(base + "testDepth.tiff"));
// sqrt
ortho.endIntegration();
ortho.saveDepth(qPrintable("final_depth.tif"));


//depthCam.camera.loadXml(orientationXmlPath);
//depthCam.loadDepth(qPrintable(tiffFile.absoluteFilePath()));

//ortho.saveDepth(qPrintable(base + "testDepth.tiff"));
//ortho.computeGaussianWeightedGrid(qPrintable(point_txt));


//int pixelX = 165;
//int pixelY = 144;
//float pixelZ = 4.5;
ortho.projectToCameraDepthMap(depthCam.camera, outputPath);

Eigen::Matrix3f rotationMatrix;
Eigen::Vector3f center;

// depth.loadMask(qPrintable(maskPath));
// depth.saveDepth(qPrintable(depthmapPath));
// depth.saveMask(qPrintable(maskPath));
//QString maskObjPath = base + "testcenterRel_copia/photogrammetry/mask.obj";
ortho.saveObj(qPrintable(base + "testcenterRel_copia/photogrammetry/depthmap_projectL05C13.obj"));
ortho.saveObj(qPrintable(base + "depthmap_projectL05C13.obj"));


//depth.depthIntegrateNormals();
Expand Down

0 comments on commit 2fcf01f

Please sign in to comment.