Commit d5246a33 authored by Lukas Werner's avatar Lukas Werner
Browse files

Moved raytracing code into its own files

parent cf94a536
Pipeline #6859 canceled with stage
in 14 seconds
......@@ -32,8 +32,7 @@
#include <pe/vtk/SphereVtkOutput.h>
#include "vtk/VTKOutput.h"
#include <pe/raytracing/Intersects.h>
#include <pe/raytracing/Ray.h>
#include <pe/raytracing/Raytracer.h>
#include <core/mpi/all.h>
......@@ -47,232 +46,6 @@ using namespace walberla::pe::raytracing;
typedef boost::tuple<Sphere, Plane, Box> BodyTypeTuple ;
//! [BodyTypeTuple]
struct BodyIntersectionInfo {
size_t imageX; // viewing plane pixel coordinates to ...
size_t imageY; // ... identify ray by pixel it intersected
walberla::id_t bodySystemID; // body which was hit
real_t t; // distance from camera to intersection point on body
};
struct Coordinates {
size_t x;
size_t y;
};
struct CoordinatesComparator {
bool operator() (const Coordinates& lhs, const Coordinates& rhs) const {
// standard lexicographical ordering
return (lhs.x < rhs.x) || (lhs.x == rhs.x && lhs.y < rhs.y);
}
};
real_t deg2rad(real_t deg) {
return deg * math::M_PI / real_t(180.0);
}
void writeTBufferToFile(const std::map<Coordinates, real_t, CoordinatesComparator>& buffer, const std::map<Coordinates, walberla::id_t, CoordinatesComparator> idBuffer, const size_t width, const size_t height, const std::string& fileName) {
real_t t_max = 1;
real_t t_min = INFINITY;
walberla::id_t maxId = 0;
for (size_t x = 0; x < width; x++) {
for (size_t y = 0; y < height; y++) {
Coordinates c = {x, y};
real_t t = buffer.at(c);
if (t > t_max && !realIsIdentical(t, INFINITY)) {
t_max = t;
}
if (t < t_min) {
t_min = t;
}
if (idBuffer.at(c) > maxId) {
maxId = idBuffer.at(c);
}
}
}
if (realIsIdentical(t_min, INFINITY)) t_min = 0;
//std::map<walberla::id_t, char> idToColors;
std::ofstream ofs(fileName, std::ios::out | std::ios::binary);
ofs << "P6\n" << width << " " << height << "\n255\n";
for (size_t y = height-1; y > 0; y--) {
for (size_t x = 0; x < width; x++) {
Coordinates c = {x, y};
char r = 0, g = 0, b = 0;
real_t t = buffer.at(c);
//walberla::id_t id = idBuffer[i];
if (realIsIdentical(t, INFINITY)) {
r = g = b = (char)255;
} else {
//r = g = b = (char)(255 * (t-t_min)/t_max);
//b = (char)(255 * (real_t(id)/real_t(maxId)));
//if (b > (char)250) b = (char)250;
//if (idToColors.count(id) == 0) {
// idToColors.insert(std::make_pair(id, math::intRandom(0, 240)));
//}
r = g = b = (char)(200 * ((t-t_min)/(t_max-t_min)));
//b = (char)(200 * ((t-t_min)/t_max));
//r = (char)(200 * ((t-t_min)/t_max));
}
ofs << r << g << b;
}
}
ofs.close();
}
//forest speichert alle blöcke meines prozesses
void rayTrace (shared_ptr<BlockForest> forest, BlockDataID storageID) {
// - settings
size_t pixelsHorizontal = 640;
size_t pixelsVertical = 480;
real_t fov_vertical = 49.13; // in degrees, in vertical direction
// camera settings
Vec3 cameraPosition(-25,10,10); // -5,0,0 for testing, -25,10,10 for simulation
Vec3 lookAtPoint(-5,10,10); // 1,0,0 for testing, -5,10,10 for simulation
Vec3 upVector(0,0,1);
// - viewing plane construction
// eye cos setup
Vec3 n = (cameraPosition - lookAtPoint).getNormalized(); // normal vector of viewing plane
Vec3 u = (upVector % n).getNormalized(); // u and ...
Vec3 v = n % u; // ... v span the viewing plane
// image plane setup
real_t d = (cameraPosition - lookAtPoint).length(); // distance from camera to viewing plane
real_t aspectRatio = real_t(pixelsHorizontal) / real_t(pixelsVertical);
real_t imagePlaneHeight = tan(deg2rad(fov_vertical)/real_t(2.)) * real_t(2.) * d;
real_t imagePlaneWidth = imagePlaneHeight * aspectRatio;
Vec3 imagePlaneOrigin = lookAtPoint - u*imagePlaneWidth/real_t(2.) - v*imagePlaneHeight/real_t(2.);
real_t pixelWidth = imagePlaneWidth / real_c(pixelsHorizontal);
real_t pixelHeight = imagePlaneHeight / real_c(pixelsVertical);
// - raytracing
std::map<Coordinates, real_t, CoordinatesComparator> tBuffer; // t values for each pixel
std::map<Coordinates, walberla::id_t, CoordinatesComparator> idBuffer; // ids of the intersected body for each pixel
std::vector<BodyIntersectionInfo> intersections; // contains for each pixel information about an intersection, if existent
std::map<Coordinates, BodyIntersectionInfo, CoordinatesComparator> localPixelIntersectionMap;
real_t t, t_closest;
walberla::id_t id_closest;
RigidBody* body_closest = NULL;
Ray ray(cameraPosition, Vec3(1,0,0));
IntersectsFunctor func(ray, t);
for (size_t x = 0; x < pixelsHorizontal; x++) {
for (size_t y = 0; y < pixelsVertical; y++) {
//WALBERLA_LOG_INFO(x << "/" << y);
Vec3 pixelLocation = imagePlaneOrigin + u*(real_c(x)+real_t(0.5))*pixelWidth + v*(real_c(y)+real_t(0.5))*pixelHeight;
Vec3 direction = (pixelLocation - cameraPosition).getNormalized();
ray.setDirection(direction);
t_closest = INFINITY;
id_closest = 0;
body_closest = NULL;
for (auto blockIt = forest->begin(); blockIt != forest->end(); ++blockIt) {
// blockIt->getAABB();
/*const AABB& blockAabb = blockIt->getAABB();
if (!intersects(blockAabb, ray, t)) {
continue;
}*/
for (auto bodyIt = LocalBodyIterator::begin(*blockIt, storageID); bodyIt != LocalBodyIterator::end(); ++bodyIt) {
bool intersects = SingleCast<BodyTypeTuple, IntersectsFunctor, bool>::execute(*bodyIt, func);
if (intersects && t < t_closest) {
// body was shot by ray and currently closest to camera
t_closest = t;
id_closest = bodyIt->getID();
body_closest = *bodyIt;
}
}
}
//std::cout << (t_closest != INFINITY ? size_t(t_closest) : 0) << " ";
Coordinates c = {
x,
y
};
if (!realIsIdentical(t_closest, INFINITY) && body_closest != NULL) {
BodyIntersectionInfo intersectionInfo = {
x,
y,
body_closest->getSystemID(),
t_closest
};
intersections.push_back(intersectionInfo);
localPixelIntersectionMap[c] = intersectionInfo;
}
tBuffer[c] = t_closest;
idBuffer[c] = id_closest;
}
//std::cout << std::endl;
}
// intersections synchronisieren
mpi::SendBuffer sendBuffer;
for (auto& info: intersections) {
sendBuffer << info.imageX << info.imageY << info.bodySystemID << info.t;
}
std::vector<BodyIntersectionInfo> gatheredIntersections;
std::map<walberla::id_t, bool> visibleBodyIDs;
//std::map<Coordinates, BodyIntersectionInfo, CoordinatesComparator> pixelIntersectionMap;
mpi::RecvBuffer recvBuffer;
mpi::allGathervBuffer(sendBuffer, recvBuffer);
while (!recvBuffer.isEmpty()) {
BodyIntersectionInfo info;
recvBuffer >> info.imageX;
recvBuffer >> info.imageY;
recvBuffer >> info.bodySystemID;
recvBuffer >> info.t;
Coordinates c = {
info.imageX,
info.imageY
};
/*if (pixelIntersectionMap.find(c) == pixelIntersectionMap.end()) {
// map didnt contain coordinates
pixelIntersectionMap.insert(std::make_pair(c, info));
} else {
// map already contains info at coordinates, check if current info is closer
BodyIntersectionInfo& existingInfo = pixelIntersectionMap.at(c);
if (existingInfo.t < info.t) {
pixelIntersectionMap[c] = info;
}
}*/
auto it = localPixelIntersectionMap.find(c);
if (it != localPixelIntersectionMap.end()) {
// there was a local hit at coordinate c
BodyIntersectionInfo& localInfo = localPixelIntersectionMap.at(c);
if (localInfo.t < info.t) {
localPixelIntersectionMap.erase(it);
}
}
//gatheredIntersections.push_back(info);
}
for (auto& info: localPixelIntersectionMap) {
visibleBodyIDs[info.second.bodySystemID] = true;
}
WALBERLA_LOG_INFO("#particles visible: " << visibleBodyIDs.size());
//WALBERLA_LOG_INFO_ON_ROOT("#gatheredIntersections: " << gatheredIntersections.size());
#ifdef __APPLE__
mpi::MPIRank rank = mpi::MPIManager::instance()->rank();
writeTBufferToFile(tBuffer, idBuffer, pixelsHorizontal, pixelsVertical, "/Users/ng/Desktop/walberla/tbuffer_" + std::to_string(rank) + ".ppm");
#endif
}
void testRayTracing () {
shared_ptr<BodyStorage> globalBodyStorage = make_shared<BodyStorage>();
shared_ptr<BlockForest> forest = createBlockForest(AABB(0,-5,-5,10,5,5), Vec3(1,1,1), Vec3(false, false, false));
......@@ -286,7 +59,8 @@ void testRayTracing () {
createBox(*globalBodyStorage, *forest, storageID, 7, Vec3(5,-4,3), Vec3(2,2,2));
//createSphere(*globalBodyStorage, *forest, storageID, 5, Vec3(1,0,0), real_t(0.1));
rayTrace(forest, storageID);
//Raytracer raytracer(forest, storageID, uint8_t(640), uint8_t(480), 49.13, Vec3(-5,0,0), Vec3(-1,0,0), Vec3(0,0,1));
//raytracer.rayTrace<BodyTypeTuple>(0);
}
int main( int argc, char ** argv )
......@@ -356,6 +130,9 @@ int main( int argc, char ** argv )
auto ccdID = forest->addBlockData(ccd::createHashGridsDataHandling( globalBodyStorage, storageID ), "CCD");
auto fcdID = forest->addBlockData(fcd::createGenericFCDDataHandling<BodyTypeTuple, fcd::AnalyticCollideFunctor>(), "FCD");
//! [AdditionalBlockData]
WALBERLA_LOG_INFO_ON_ROOT("*** RAYTRACER ***");
Raytracer raytracer(forest, storageID, uint8_t(640), uint8_t(480), 49.13, Vec3(-25,10,10), Vec3(-1,10,10), Vec3(0,0,1));
WALBERLA_LOG_INFO_ON_ROOT("*** INTEGRATOR ***");
//! [Integrator]
......@@ -439,7 +216,7 @@ int main( int argc, char ** argv )
//! [PostProcessing]
WALBERLA_LOG_INFO_ON_ROOT("*** RAYTRACING - START ***");
rayTrace(forest, storageID);
raytracer.rayTrace<BodyTypeTuple>(0);
WALBERLA_LOG_INFO_ON_ROOT("*** RAYTRACING - END ***");
// für einzelne sphere vtks: -> SphereVtkOutput.cpp
......
//======================================================================================================================
//
// This file is part of waLBerla. waLBerla is free software: you can
// redistribute it and/or modify it under the terms of the GNU General Public
// License as published by the Free Software Foundation, either version 3 of
// the License, or (at your option) any later version.
//
// waLBerla is distributed in the hope that it will be useful, but WITHOUT
// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
// for more details.
//
// You should have received a copy of the GNU General Public License along
// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>.
//
//! \file Raytracer.cpp
//! \author Lukas Werner
//
//======================================================================================================================
#include "Raytracer.h"
#include "Ray.h"
#include "Intersects.h"
using namespace walberla;
struct BodyIntersectionInfo {
size_t imageX; // viewing plane pixel coordinates to ...
size_t imageY; // ... identify ray by pixel it intersected
walberla::id_t bodySystemID; // body which was hit
real_t t; // distance from camera to intersection point on body
};
struct Coordinates {
size_t x;
size_t y;
};
struct CoordinatesComparator {
bool operator() (const Coordinates& lhs, const Coordinates& rhs) const {
// standard lexicographical ordering
return (lhs.x < rhs.x) || (lhs.x == rhs.x && lhs.y < rhs.y);
}
};
real_t deg2rad(real_t deg) {
return deg * math::M_PI / real_t(180.0);
}
namespace walberla {
namespace pe {
namespace raytracing {
/*!\brief Instantiation constructor for the Box class.
*
* \param forest BlockForest the raytracer operates on.
* \param storageID Storage ID of the block data storage the raytracer operates on.
* \param pixelsHorizontal Horizontal amount of pixels of the generated image.
* \param pixelsVertical Vertical amount of pixels of the generated image.
* \param fov_vertical Vertical field-of-view of the camera.
* \param cameraPosition Position of the camera in the global world frame.
* \param lookAtPoint Point the camera looks at in the global world frame.
* \param upVector Vector indicating the upwards direction of the camera.
*/
Raytracer::Raytracer(const shared_ptr<BlockStorage>& forest, BlockDataID storageID,
uint8_t pixelsHorizontal, uint8_t pixelsVertical,
real_t fov_vertical,
const Vec3& cameraPosition, const Vec3& lookAtPoint, const Vec3& upVector)
: forest_(forest), storageID_(storageID),
pixelsHorizontal_(pixelsHorizontal), pixelsVertical_(pixelsVertical),
fov_vertical_(fov_vertical),
cameraPosition_(cameraPosition), lookAtPoint_(lookAtPoint), upVector_(upVector)
{
// eye coordinate system setup
n = (cameraPosition_ - lookAtPoint_).getNormalized();
u = (upVector_ % n).getNormalized();
v = n % u;
// viewing plane setup
d = (cameraPosition_ - lookAtPoint_).length();
aspectRatio = real_t(pixelsHorizontal_) / real_t(pixelsVertical_);
viewingPlaneHeight = tan(deg2rad(fov_vertical_)/real_t(2.)) * real_t(2.) * d;
viewingPlaneWidth = viewingPlaneHeight * aspectRatio;
viewingPlaneOrigin = lookAtPoint_ - u*viewingPlaneWidth/real_t(2.) - v*viewingPlaneHeight/real_t(2.);
pixelWidth = viewingPlaneWidth / real_c(pixelsHorizontal_);
pixelHeight = viewingPlaneHeight / real_c(pixelsVertical);
}
/*!\brief Does one raytracing step.
*
* \param timestep The timestep after which the raytracing starts.
*/
template <typename BodyTypeTuple>
void Raytracer::rayTrace(const size_t timestep) const {
std::map<Coordinates, real_t, CoordinatesComparator> tBuffer; // t values for each pixel
std::map<Coordinates, walberla::id_t, CoordinatesComparator> idBuffer; // ids of the intersected body for each pixel
std::vector<BodyIntersectionInfo> intersections; // contains for each pixel information about an intersection, if existent
std::map<Coordinates, BodyIntersectionInfo, CoordinatesComparator> localPixelIntersectionMap;
real_t t, t_closest;
walberla::id_t id_closest;
RigidBody* body_closest = NULL;
Ray ray(cameraPosition_, Vec3(1,0,0));
IntersectsFunctor func(ray, t);
for (size_t x = 0; x < pixelsHorizontal_; x++) {
for (size_t y = 0; y < pixelsVertical_; y++) {
//WALBERLA_LOG_INFO(x << "/" << y);
Vec3 pixelLocation = viewingPlaneOrigin + u*(real_c(x)+real_t(0.5))*pixelWidth + v*(real_c(y)+real_t(0.5))*pixelHeight;
Vec3 direction = (pixelLocation - cameraPosition_).getNormalized();
ray.setDirection(direction);
t_closest = INFINITY;
id_closest = 0;
body_closest = NULL;
for (auto blockIt = forest_->begin(); blockIt != forest_->end(); ++blockIt) {
// blockIt->getAABB();
/*const AABB& blockAabb = blockIt->getAABB();
if (!intersects(blockAabb, ray, t)) {
continue;
}*/
for (auto bodyIt = LocalBodyIterator::begin(*blockIt, storageID_); bodyIt != LocalBodyIterator::end(); ++bodyIt) {
bool intersects = SingleCast<BodyTypeTuple, IntersectsFunctor, bool>::execute(*bodyIt, func);
if (intersects && t < t_closest) {
// body was shot by ray and currently closest to camera
t_closest = t;
id_closest = bodyIt->getID();
body_closest = *bodyIt;
}
}
}
//std::cout << (t_closest != INFINITY ? size_t(t_closest) : 0) << " ";
Coordinates c = {
x,
y
};
if (!realIsIdentical(t_closest, INFINITY) && body_closest != NULL) {
BodyIntersectionInfo intersectionInfo = {
x,
y,
body_closest->getSystemID(),
t_closest
};
intersections.push_back(intersectionInfo);
localPixelIntersectionMap[c] = intersectionInfo;
}
tBuffer[c] = t_closest;
idBuffer[c] = id_closest;
}
//std::cout << std::endl;
}
// intersections synchronisieren
mpi::SendBuffer sendBuffer;
for (auto& info: intersections) {
sendBuffer << info.imageX << info.imageY << info.bodySystemID << info.t;
}
std::vector<BodyIntersectionInfo> gatheredIntersections;
std::map<walberla::id_t, bool> visibleBodyIDs;
//std::map<Coordinates, BodyIntersectionInfo, CoordinatesComparator> pixelIntersectionMap;
mpi::RecvBuffer recvBuffer;
mpi::allGathervBuffer(sendBuffer, recvBuffer);
while (!recvBuffer.isEmpty()) {
BodyIntersectionInfo info;
recvBuffer >> info.imageX;
recvBuffer >> info.imageY;
recvBuffer >> info.bodySystemID;
recvBuffer >> info.t;
Coordinates c = {
info.imageX,
info.imageY
};
/*if (pixelIntersectionMap.find(c) == pixelIntersectionMap.end()) {
// map didnt contain coordinates
pixelIntersectionMap.insert(std::make_pair(c, info));
} else {
// map already contains info at coordinates, check if current info is closer
BodyIntersectionInfo& existingInfo = pixelIntersectionMap.at(c);
if (existingInfo.t < info.t) {
pixelIntersectionMap[c] = info;
}
}*/
auto it = localPixelIntersectionMap.find(c);
if (it != localPixelIntersectionMap.end()) {
// there was a local hit at coordinate c
BodyIntersectionInfo& localInfo = localPixelIntersectionMap.at(c);
if (localInfo.t < info.t) {
localPixelIntersectionMap.erase(it);
}
}
//gatheredIntersections.push_back(info);
}
for (auto& info: localPixelIntersectionMap) {
visibleBodyIDs[info.second.bodySystemID] = true;
}
WALBERLA_LOG_INFO("#particles visible: " << visibleBodyIDs.size());
//WALBERLA_LOG_INFO_ON_ROOT("#gatheredIntersections: " << gatheredIntersections.size());
#ifdef __APPLE__
//mpi::MPIRank rank = mpi::MPIManager::instance()->rank();
//writeTBufferToFile(tBuffer, idBuffer, pixelsHorizontal_, pixelsVertical_, "/Users/ng/Desktop/walberla/tbuffer_" + std::to_string(rank) + ".ppm");
#endif
}
}
}
}
//======================================================================================================================
//
// This file is part of waLBerla. waLBerla is free software: you can
// redistribute it and/or modify it under the terms of the GNU General Public
// License as published by the Free Software Foundation, either version 3 of
// the License, or (at your option) any later version.
//
// waLBerla is distributed in the hope that it will be useful, but WITHOUT
// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
// for more details.
//
// You should have received a copy of the GNU General Public License along
// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>.
//
//! \file Raytracer.h
//! \author Lukas Werner
//
//======================================================================================================================
#pragma once
#include <pe/basic.h>
#include <pe/Types.h>
#include <core/math/Vector3.h>
#include <core/mpi/all.h>
namespace walberla {
namespace pe {
namespace raytracing {
class Raytracer {
public:
/*!\name Constructors */
//@{
explicit Raytracer(const shared_ptr<BlockStorage> & forest, BlockDataID storageID,
uint8_t pixelsHorizontal, uint8_t pixelsVertical,
real_t fov_vertical,
const Vec3& cameraPosition, const Vec3& lookAtPoint, const Vec3& upVector);
//@}
private:
/*!\name Member variables */
//@{
const shared_ptr<BlockStorage>& forest_; //!< The BlockForest the raytracer operates on.
BlockDataID storageID_; /*!< The storage ID of the block data storage the raytracer operates
on.*/
uint8_t pixelsHorizontal_; //!< The horizontal amount of pixels of the generated image.
uint8_t pixelsVertical_; //!< The vertical amount of pixels of the generated image.
real_t fov_vertical_; //!< The vertical field-of-view of the camera.
Vec3 cameraPosition_; //!< The position of the camera in the global world frame.
Vec3 lookAtPoint_; /*!< The point the camera looks at in the global world frame,
marks the center of the view plane.*/
Vec3 upVector_; //!< The vector indicating the upwards direction of the camera.
//@}
Vec3 n; // normal vector of viewing plane
Vec3 u; // u and ...
Vec3 v; // ... v span the viewing plane
real_t d; // distance from camera to viewing plane
real_t aspectRatio; // aspect ratio of the generated image and viewing plane
real_t viewingPlaneHeight; // viewing plane height
real_t viewingPlaneWidth; // viewing plane width
Vec3 viewingPlaneOrigin; // origin of the viewing plane
real_t pixelWidth; // width of a pixel of the generated image in the viewing plane
real_t pixelHeight; // height of a pixel of the generated image in the viewing plane
public:
/*!\name Get functions */
//@{
inline uint8_t getPixelsHorizontal() const;
inline uint8_t getPixelsVertical() const;
inline real_t getFOVVertical() const;
inline const Vec3& getCameraPosition() const;
inline const Vec3& getLookAtPoint() const;
inline const Vec3& getUpVector() const;
//@}
/*!\name Functions */
//@{
template <typename BodyTypeTuple>
void rayTrace(const size_t timestep) const;
//@}
};
/*!\brief Returns the horizontal amount of pixels of the generated image.
*
* \return The horizontal amount of pixels of the generated image.
*/
inline uint8_t Raytracer::getPixelsHorizontal() const {
return pixelsHorizontal_;
}
/*!\brief Returns the vertical amount of pixels of the generated image.
*
* \return The vertical amount of pixels of the generated image.
*/
inline uint8_t Raytracer::getPixelsVertical() const {
return pixelsVertical_;
}
/*!\brief Returns the vertical field-of-view of the camera.
*
* \return The vertical field-of-view of the camera.
*/
inline real_t Raytracer::getFOVVertical() const {
return fov_vertical_;
}
/*!\brief Returns the position of the camera in the global world frame.
*
* \return The position of the camera.