From 481375a6b766ffd7a7738e48e2402b6ee3fb501a Mon Sep 17 00:00:00 2001 From: robin Date: Mon, 6 May 2024 11:11:12 +0200 Subject: [PATCH] initial rough working --- advection/src/CMakeLists.txt | 42 ----------- advection/src/UVGrid.cpp | 66 ----------------- vtk/src/CMakeLists.txt | 16 +++++ .../src/advection}/AdvectionKernel.h | 2 +- .../src/advection}/EulerAdvectionKernel.cpp | 0 .../src/advection}/EulerAdvectionKernel.h | 0 .../src/advection}/RK4AdvectionKernel.cpp | 0 .../src/advection}/RK4AdvectionKernel.h | 0 vtk/src/advection/UVGrid.cpp | 66 +++++++++++++++++ {advection/src => vtk/src/advection}/UVGrid.h | 0 {advection/src => vtk/src/advection}/Vel.cpp | 0 {advection/src => vtk/src/advection}/Vel.h | 0 .../src => vtk/src/advection}/interpolate.cpp | 0 .../src => vtk/src/advection}/interpolate.h | 0 {advection/src => vtk/src/advection}/main.cpp | 0 .../src => vtk/src/advection}/readdata.cpp | 0 .../src => vtk/src/advection}/readdata.h | 0 vtk/src/layers/EGlyphLayer.cpp | 71 ++++++++----------- vtk/src/layers/EGlyphLayer.h | 5 +- vtk/src/layers/LGlyphLayer.cpp | 6 +- vtk/src/layers/LGlyphLayer.h | 4 +- vtk/src/main.cpp | 12 +++- 22 files changed, 131 insertions(+), 159 deletions(-) delete mode 100644 advection/src/CMakeLists.txt delete mode 100644 advection/src/UVGrid.cpp rename {advection/src => vtk/src/advection}/AdvectionKernel.h (94%) rename {advection/src => vtk/src/advection}/EulerAdvectionKernel.cpp (100%) rename {advection/src => vtk/src/advection}/EulerAdvectionKernel.h (100%) rename {advection/src => vtk/src/advection}/RK4AdvectionKernel.cpp (100%) rename {advection/src => vtk/src/advection}/RK4AdvectionKernel.h (100%) create mode 100644 vtk/src/advection/UVGrid.cpp rename {advection/src => vtk/src/advection}/UVGrid.h (100%) rename {advection/src => vtk/src/advection}/Vel.cpp (100%) rename {advection/src => vtk/src/advection}/Vel.h (100%) rename {advection/src => vtk/src/advection}/interpolate.cpp (100%) rename {advection/src => vtk/src/advection}/interpolate.h (100%) rename {advection/src => vtk/src/advection}/main.cpp (100%) rename {advection/src => vtk/src/advection}/readdata.cpp (100%) rename {advection/src => vtk/src/advection}/readdata.h (100%) diff --git a/advection/src/CMakeLists.txt b/advection/src/CMakeLists.txt deleted file mode 100644 index b8c8758..0000000 --- a/advection/src/CMakeLists.txt +++ /dev/null @@ -1,42 +0,0 @@ -cmake_minimum_required (VERSION 3.28) -project (Advection) - -set(CMAKE_CXX_STANDARD 23) -set(CMAKE_CXX_STANDARD_REQUIRED ON) - -set(CMAKE_EXPORT_COMPILE_COMMANDS ON) - -find_package(netCDF REQUIRED) - -add_executable(Advection main.cpp - readdata.cpp - readdata.h - interpolate.cpp - interpolate.h - UVGrid.cpp - UVGrid.h - Vel.h - Vel.cpp - AdvectionKernel.h - EulerAdvectionKernel.cpp - EulerAdvectionKernel.h - RK4AdvectionKernel.cpp - RK4AdvectionKernel.h -) - -execute_process( - COMMAND nc-config --includedir - OUTPUT_VARIABLE NETCDF_INCLUDE_DIR - OUTPUT_STRIP_TRAILING_WHITESPACE -) - -execute_process( - COMMAND ncxx4-config --libdir - OUTPUT_VARIABLE NETCDFCXX_LIB_DIR - OUTPUT_STRIP_TRAILING_WHITESPACE -) - -target_include_directories(Advection PUBLIC ${netCDF_INCLUDE_DIR}) - -find_library(NETCDF_LIB NAMES netcdf-cxx4 netcdf_c++4 PATHS ${NETCDFCXX_LIB_DIR} NO_DEFAULT_PATH) -target_link_libraries(Advection ${NETCDF_LIB}) diff --git a/advection/src/UVGrid.cpp b/advection/src/UVGrid.cpp deleted file mode 100644 index 1febbfe..0000000 --- a/advection/src/UVGrid.cpp +++ /dev/null @@ -1,66 +0,0 @@ -#include - -#include "UVGrid.h" -#include "readdata.h" - -#define sizeError2 "The sizes of the hydrodynamic data files are different" -#define sizeError "The sizes of the hydrodynamicU or -V files does not correspond with the sizes of the grid file" - -using namespace std; - -UVGrid::UVGrid() { - auto us = readHydrodynamicU(); - auto vs = readHydrodynamicV(); - if (us.size() != vs.size()) { - throw domain_error(sizeError2); - } - - tie(times, lats, lons) = readGrid(); - - timeSize = times.size(); - latSize = lats.size(); - lonSize = lons.size(); - - size_t gridSize = timeSize * latSize * lonSize; - if (gridSize != us.size()) { - throw domain_error(sizeError); - } - - uvData.reserve(gridSize); - - for (auto vel: views::zip(us, vs)) { - uvData.push_back(Vel(vel)); - } -} - -const Vel &UVGrid::operator[](size_t timeIndex, size_t latIndex, size_t lonIndex) const { - if(timeIndex < 0 or timeIndex >= timeSize - or latIndex < 0 or latIndex >= latSize - or lonIndex < 0 or lonIndex >= lonSize) { - throw std::out_of_range("Index out of bounds"); - } - size_t index = timeIndex * (latSize * lonSize) + latIndex * lonSize + lonIndex; - return uvData[index]; -} - -double UVGrid::lonStep() const { - return lons[1] - lons[0]; -} - -double UVGrid::latStep() const { - return lats[1] - lats[0]; -} - -int UVGrid::timeStep() const { - return times[1] - times[0]; -} - -void UVGrid::streamSlice(ostream &os, size_t t) { - for (int x = 0; x < latSize; x++) { - for (int y = 0; y < lonSize; y++) { - auto vel = (*this)[t,x,y]; - os << vel << " "; - } - os << endl; - } -} \ No newline at end of file diff --git a/vtk/src/CMakeLists.txt b/vtk/src/CMakeLists.txt index ef2a3aa..cf0b9d3 100644 --- a/vtk/src/CMakeLists.txt +++ b/vtk/src/CMakeLists.txt @@ -2,6 +2,9 @@ cmake_minimum_required(VERSION 3.12 FATAL_ERROR) project(VtkBase) +set(CMAKE_CXX_STANDARD 23) +set(CMAKE_CXX_STANDARD_REQUIRED ON) + set(CMAKE_EXPORT_COMPILE_COMMANDS ON) find_package(VTK COMPONENTS @@ -47,6 +50,19 @@ add_executable(VtkBase MACOSX_BUNDLE main.cpp commands/SpawnPointCallback.h commands/SpawnPointCallback.cpp CartographicTransformation.cpp + advection/AdvectionKernel.h + advection/EulerAdvectionKernel.cpp + advection/EulerAdvectionKernel.h + advection/interpolate.cpp + advection/interpolate.h + advection/readdata.cpp + advection/readdata.h + advection/RK4AdvectionKernel.cpp + advection/RK4AdvectionKernel.h + advection/UVGrid.cpp + advection/UVGrid.h + advection/Vel.cpp + advection/Vel.h ) execute_process( diff --git a/advection/src/AdvectionKernel.h b/vtk/src/advection/AdvectionKernel.h similarity index 94% rename from advection/src/AdvectionKernel.h rename to vtk/src/advection/AdvectionKernel.h index 0170c79..d8d7674 100644 --- a/advection/src/AdvectionKernel.h +++ b/vtk/src/advection/AdvectionKernel.h @@ -11,7 +11,7 @@ */ class AdvectionKernel { public: - const static int DT = 60 * 15; // 60 sec/min * 15 mins + const static int DT = 15 * 60 * 15; // 60 sec/min * 15 mins /** * This function must take a time, latitude and longitude of a particle and must output * a new latitude and longitude after being advected once for AdvectionKernel::DT time as defined above. diff --git a/advection/src/EulerAdvectionKernel.cpp b/vtk/src/advection/EulerAdvectionKernel.cpp similarity index 100% rename from advection/src/EulerAdvectionKernel.cpp rename to vtk/src/advection/EulerAdvectionKernel.cpp diff --git a/advection/src/EulerAdvectionKernel.h b/vtk/src/advection/EulerAdvectionKernel.h similarity index 100% rename from advection/src/EulerAdvectionKernel.h rename to vtk/src/advection/EulerAdvectionKernel.h diff --git a/advection/src/RK4AdvectionKernel.cpp b/vtk/src/advection/RK4AdvectionKernel.cpp similarity index 100% rename from advection/src/RK4AdvectionKernel.cpp rename to vtk/src/advection/RK4AdvectionKernel.cpp diff --git a/advection/src/RK4AdvectionKernel.h b/vtk/src/advection/RK4AdvectionKernel.h similarity index 100% rename from advection/src/RK4AdvectionKernel.h rename to vtk/src/advection/RK4AdvectionKernel.h diff --git a/vtk/src/advection/UVGrid.cpp b/vtk/src/advection/UVGrid.cpp new file mode 100644 index 0000000..f35aa72 --- /dev/null +++ b/vtk/src/advection/UVGrid.cpp @@ -0,0 +1,66 @@ +#include + +#include "UVGrid.h" +#include "readdata.h" + +#define sizeError2 "The sizes of the hydrodynamic data files are different" +#define sizeError "The sizes of the hydrodynamicU or -V files does not correspond with the sizes of the grid file" + +using namespace std; + +UVGrid::UVGrid() { + auto us = readHydrodynamicU(); + auto vs = readHydrodynamicV(); + if (us.size() != vs.size()) { + throw domain_error(sizeError2); + } + + tie(times, lats, lons) = readGrid(); + + timeSize = times.size(); + latSize = lats.size(); + lonSize = lons.size(); + + size_t gridSize = timeSize * latSize * lonSize; + if (gridSize != us.size()) { + throw domain_error(sizeError); + } + + uvData.reserve(gridSize); + + for (auto vel: views::zip(us, vs)) { + uvData.push_back(Vel(vel)); + } +} + +const Vel &UVGrid::operator[](size_t timeIndex, size_t latIndex, size_t lonIndex) const { + if (timeIndex < 0 or timeIndex >= timeSize + or latIndex < 0 or latIndex >= latSize + or lonIndex < 0 or lonIndex >= lonSize) { + throw std::out_of_range("Index out of bounds"); + } + size_t index = timeIndex * (latSize * lonSize) + latIndex * lonSize + lonIndex; + return uvData[index]; +} + +double UVGrid::lonStep() const { + return lons[1] - lons[0]; +} + +double UVGrid::latStep() const { + return lats[1] - lats[0]; +} + +int UVGrid::timeStep() const { + return times[1] - times[0]; +} + +void UVGrid::streamSlice(ostream &os, size_t t) { + for (int x = 0; x < latSize; x++) { + for (int y = 0; y < lonSize; y++) { + auto vel = (*this)[t, x, y]; + os << vel << " "; + } + os << endl; + } +} \ No newline at end of file diff --git a/advection/src/UVGrid.h b/vtk/src/advection/UVGrid.h similarity index 100% rename from advection/src/UVGrid.h rename to vtk/src/advection/UVGrid.h diff --git a/advection/src/Vel.cpp b/vtk/src/advection/Vel.cpp similarity index 100% rename from advection/src/Vel.cpp rename to vtk/src/advection/Vel.cpp diff --git a/advection/src/Vel.h b/vtk/src/advection/Vel.h similarity index 100% rename from advection/src/Vel.h rename to vtk/src/advection/Vel.h diff --git a/advection/src/interpolate.cpp b/vtk/src/advection/interpolate.cpp similarity index 100% rename from advection/src/interpolate.cpp rename to vtk/src/advection/interpolate.cpp diff --git a/advection/src/interpolate.h b/vtk/src/advection/interpolate.h similarity index 100% rename from advection/src/interpolate.h rename to vtk/src/advection/interpolate.h diff --git a/advection/src/main.cpp b/vtk/src/advection/main.cpp similarity index 100% rename from advection/src/main.cpp rename to vtk/src/advection/main.cpp diff --git a/advection/src/readdata.cpp b/vtk/src/advection/readdata.cpp similarity index 100% rename from advection/src/readdata.cpp rename to vtk/src/advection/readdata.cpp diff --git a/advection/src/readdata.h b/vtk/src/advection/readdata.h similarity index 100% rename from advection/src/readdata.h rename to vtk/src/advection/readdata.h diff --git a/vtk/src/layers/EGlyphLayer.cpp b/vtk/src/layers/EGlyphLayer.cpp index cb91825..afcb4b3 100644 --- a/vtk/src/layers/EGlyphLayer.cpp +++ b/vtk/src/layers/EGlyphLayer.cpp @@ -11,43 +11,20 @@ #include #include #include -#include #include #include "../CartographicTransformation.h" +#include "../advection/readdata.h" +#include "../advection/interpolate.h" -using namespace netCDF; using namespace std; -template -vector getVarVector(const NcVar &var) { - int length = 1; - for (NcDim dim : var.getDims()) { - length *= dim.getSize(); - } - - vector vec(length); - - var.getVar(vec.data()); - - return vec; -} - -tuple, vector, vector> readGrid() { - netCDF::NcFile data("../../../../data/grid.h5", netCDF::NcFile::read); - multimap< string, NcVar > vars = data.getVars(); - vector time = getVarVector(vars.find("times")->second); - vector longitude = getVarVector(vars.find("longitude")->second); - vector latitude = getVarVector(vars.find("latitude")->second); - - return {time, latitude, longitude}; -} - - -EGlyphLayer::EGlyphLayer() { +EGlyphLayer::EGlyphLayer(std::shared_ptr uvGrid) { this->ren = vtkSmartPointer::New(); this->ren->SetLayer(1); this->ren->InteractiveOff(); + this->uvGrid = uvGrid; + this->data = vtkSmartPointer::New(); this->direction = vtkSmartPointer::New(); this->direction->SetName("direction"); @@ -57,25 +34,28 @@ EGlyphLayer::EGlyphLayer() { void EGlyphLayer::readCoordinates() { vtkNew points; - auto [times, lats, lons] = readGrid(); // FIXME: import Robin's readData function and use it - this->numLats = lats.size(); - this->numLons = lons.size(); + this->numLats = uvGrid->lats.size(); + this->numLons = uvGrid->lons.size(); this->direction->SetNumberOfComponents(3); - this->direction->SetNumberOfTuples(numLats*numLons); - points->Allocate(numLats*numLons); + this->direction->SetNumberOfTuples(numLats * numLons); + points->Allocate(numLats * numLons); auto camera = createNormalisedCamera(); ren->SetActiveCamera(camera); int i = 0; - for (double lat : lats) { - for (double lon : lons) { - cout << "lon: " << lon << " lat: " << lat << endl; - direction->SetTuple3(i, 0.45, 0.90, 0); //FIXME: read this info from file + int latIndex = 0; + for (double lat: uvGrid->lats) { + int lonIndex = 0; + for (double lon: uvGrid->lons) { + auto [u, v] = (*uvGrid)[0, latIndex, lonIndex]; + direction->SetTuple3(i, 2*v, 2*u, 0); points->InsertPoint(i++, lon, lat, 0); // see also https://vtk.org/doc/nightly/html/classvtkPolyDataMapper2D.html + lonIndex++; } + latIndex++; } this->data->SetPoints(points); this->data->GetPointData()->AddArray(this->direction); @@ -94,8 +74,8 @@ void EGlyphLayer::readCoordinates() { glyph2D->SetInputConnection(transformFilter->GetOutputPort()); glyph2D->OrientOn(); glyph2D->ClampingOn(); - glyph2D->SetScaleModeToScaleByVector(); - glyph2D->SetVectorModeToUseVector(); + glyph2D->SetScaleModeToScaleByVector(); + glyph2D->SetVectorModeToUseVector(); glyph2D->Update(); // vtkNew coordinate; @@ -109,16 +89,21 @@ void EGlyphLayer::readCoordinates() { vtkNew actor; actor->SetMapper(mapper); - actor->GetProperty()->SetColor(0,0,0); + actor->GetProperty()->SetColor(0, 0, 0); actor->GetProperty()->SetOpacity(0.2); - this->ren->AddActor(actor) ; + this->ren->AddActor(actor); } void EGlyphLayer::updateData(int t) { - for (int i=0; i < numLats*numLons; i++) { - this->direction->SetTuple3(i, std::cos(t), std::sin(t), 0); // FIXME: fetch data from file. + int i = 0; + for (int lat = 0; lat < uvGrid->latSize; lat++) { + for (int lon = 0; lon < uvGrid->lonSize; lon++) { + auto [u, v] = (*uvGrid)[t/3600, lat, lon]; + this->direction->SetTuple3(i, 5*v, 5*u, 0); // FIXME: fetch data from file. + i++; + } } this->direction->Modified(); } diff --git a/vtk/src/layers/EGlyphLayer.h b/vtk/src/layers/EGlyphLayer.h index 8631f32..1dcc56a 100644 --- a/vtk/src/layers/EGlyphLayer.h +++ b/vtk/src/layers/EGlyphLayer.h @@ -4,6 +4,8 @@ #include "Layer.h" #include +#include "../advection/UVGrid.h" + /** Implements the Layer class for the case of a Eulerian visualization. * Specifically, this class models the eulerian flow-fields of the simulation using the 'glyph' mark and 'direction' and 'form' channels to denote direction and strength of velocities. */ @@ -11,6 +13,7 @@ class EGlyphLayer : public Layer { private: vtkSmartPointer data; vtkSmartPointer direction; + std::shared_ptr uvGrid; int numLats; int numLons; @@ -22,7 +25,7 @@ private: public: /** Constructor. */ - EGlyphLayer(); + EGlyphLayer(std::shared_ptr uvGrid); /** updates the glyphs to reflect the given timestamp in the dataset. * @param t : the time at which to fetch the data. diff --git a/vtk/src/layers/LGlyphLayer.cpp b/vtk/src/layers/LGlyphLayer.cpp index 845d628..9b0840a 100644 --- a/vtk/src/layers/LGlyphLayer.cpp +++ b/vtk/src/layers/LGlyphLayer.cpp @@ -31,7 +31,7 @@ vtkSmartPointer LGlyphLayer::createSpawnPointCallback() { // // TODO: modelling all this in vtkClasses is workable, but ideally i would want to work with a native C++ class. See if this is doable and feasible. -LGlyphLayer::LGlyphLayer() { +LGlyphLayer::LGlyphLayer(std::unique_ptr advectionKernel) { this->ren = vtkSmartPointer::New(); this->ren->SetLayer(2); @@ -39,6 +39,8 @@ LGlyphLayer::LGlyphLayer() { this->data = vtkSmartPointer::New(); this->data->SetPoints(this->points); + advector = std::move(advectionKernel); + auto camera = createNormalisedCamera(); ren->SetActiveCamera(camera); @@ -88,7 +90,7 @@ void LGlyphLayer::updateData(int t) { double point[3]; for (vtkIdType n = 0; n < this->points->GetNumberOfPoints(); n++) { this->points->GetPoint(n, point); - auto [xNew, yNew] = advect(n, point[0], point[1]); + auto [yNew, xNew] = advector->advect(t, point[1], point[0]); this->points->SetPoint(n, xNew, yNew, 0); } this->points->Modified(); diff --git a/vtk/src/layers/LGlyphLayer.h b/vtk/src/layers/LGlyphLayer.h index 22993b4..1a8dc0f 100644 --- a/vtk/src/layers/LGlyphLayer.h +++ b/vtk/src/layers/LGlyphLayer.h @@ -2,6 +2,7 @@ #define LGLYPHLAYER_H #include "Layer.h" +#include "../advection/AdvectionKernel.h" #include "../commands/SpawnPointCallback.h" #include #include @@ -13,13 +14,14 @@ class LGlyphLayer : public Layer { private: vtkSmartPointer points; vtkSmartPointer data; + std::unique_ptr advector; public: /** Constructor. */ - LGlyphLayer(); + LGlyphLayer(std::unique_ptr advectionKernel); /** This function spoofs a few points in the dataset. Mostly used for testing. */ diff --git a/vtk/src/main.cpp b/vtk/src/main.cpp index fb89ea4..c4ba1e7 100644 --- a/vtk/src/main.cpp +++ b/vtk/src/main.cpp @@ -7,22 +7,28 @@ #include #include #include +#include #include "layers/BackgroundImage.h" #include "layers/EGlyphLayer.h" #include "layers/LGlyphLayer.h" #include "Program.h" +#include "advection/UVGrid.h" +#include "advection/RK4AdvectionKernel.h" using namespace std; int main() { - auto l = new LGlyphLayer(); - l->spoofPoints(); + shared_ptr uvGrid = std::make_shared(); + auto kernelRK4 = make_unique(uvGrid); + + auto l = new LGlyphLayer(move(kernelRK4)); +// l->spoofPoints(); Program *program = new Program(); program->addLayer(new BackgroundImage("../../../../data/map_661-661.png")); - program->addLayer(new EGlyphLayer()); + program->addLayer(new EGlyphLayer(uvGrid)); program->addLayer(l); program->render();