diff --git a/vtk/src/advection/AdvectionKernel.h b/advection/src/AdvectionKernel.h similarity index 94% rename from vtk/src/advection/AdvectionKernel.h rename to advection/src/AdvectionKernel.h index d8d7674..0170c79 100644 --- a/vtk/src/advection/AdvectionKernel.h +++ b/advection/src/AdvectionKernel.h @@ -11,7 +11,7 @@ */ class AdvectionKernel { public: - const static int DT = 15 * 60 * 15; // 60 sec/min * 15 mins + const static int DT = 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/CMakeLists.txt b/advection/src/CMakeLists.txt new file mode 100644 index 0000000..b8c8758 --- /dev/null +++ b/advection/src/CMakeLists.txt @@ -0,0 +1,42 @@ +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/vtk/src/advection/EulerAdvectionKernel.cpp b/advection/src/EulerAdvectionKernel.cpp similarity index 100% rename from vtk/src/advection/EulerAdvectionKernel.cpp rename to advection/src/EulerAdvectionKernel.cpp diff --git a/vtk/src/advection/EulerAdvectionKernel.h b/advection/src/EulerAdvectionKernel.h similarity index 100% rename from vtk/src/advection/EulerAdvectionKernel.h rename to advection/src/EulerAdvectionKernel.h diff --git a/vtk/src/advection/RK4AdvectionKernel.cpp b/advection/src/RK4AdvectionKernel.cpp similarity index 100% rename from vtk/src/advection/RK4AdvectionKernel.cpp rename to advection/src/RK4AdvectionKernel.cpp diff --git a/vtk/src/advection/RK4AdvectionKernel.h b/advection/src/RK4AdvectionKernel.h similarity index 100% rename from vtk/src/advection/RK4AdvectionKernel.h rename to advection/src/RK4AdvectionKernel.h diff --git a/advection/src/UVGrid.cpp b/advection/src/UVGrid.cpp new file mode 100644 index 0000000..1febbfe --- /dev/null +++ b/advection/src/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/vtk/src/advection/UVGrid.h b/advection/src/UVGrid.h similarity index 100% rename from vtk/src/advection/UVGrid.h rename to advection/src/UVGrid.h diff --git a/vtk/src/advection/Vel.cpp b/advection/src/Vel.cpp similarity index 100% rename from vtk/src/advection/Vel.cpp rename to advection/src/Vel.cpp diff --git a/vtk/src/advection/Vel.h b/advection/src/Vel.h similarity index 100% rename from vtk/src/advection/Vel.h rename to advection/src/Vel.h diff --git a/vtk/src/advection/interpolate.cpp b/advection/src/interpolate.cpp similarity index 100% rename from vtk/src/advection/interpolate.cpp rename to advection/src/interpolate.cpp diff --git a/vtk/src/advection/interpolate.h b/advection/src/interpolate.h similarity index 100% rename from vtk/src/advection/interpolate.h rename to advection/src/interpolate.h diff --git a/vtk/src/advection/main.cpp b/advection/src/main.cpp similarity index 100% rename from vtk/src/advection/main.cpp rename to advection/src/main.cpp diff --git a/vtk/src/advection/readdata.cpp b/advection/src/readdata.cpp similarity index 100% rename from vtk/src/advection/readdata.cpp rename to advection/src/readdata.cpp diff --git a/vtk/src/advection/readdata.h b/advection/src/readdata.h similarity index 100% rename from vtk/src/advection/readdata.h rename to advection/src/readdata.h diff --git a/vtk/src/CMakeLists.txt b/vtk/src/CMakeLists.txt index cf0b9d3..ef2a3aa 100644 --- a/vtk/src/CMakeLists.txt +++ b/vtk/src/CMakeLists.txt @@ -2,9 +2,6 @@ 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 @@ -50,19 +47,6 @@ 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/vtk/src/advection/UVGrid.cpp b/vtk/src/advection/UVGrid.cpp deleted file mode 100644 index 3be14d9..0000000 --- a/vtk/src/advection/UVGrid.cpp +++ /dev/null @@ -1,67 +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 { - cout << "t=" << timeIndex << "latIndex=" << latIndex << "lonIndex=" << lonIndex << endl; - 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/layers/EGlyphLayer.cpp b/vtk/src/layers/EGlyphLayer.cpp index b591b54..cb91825 100644 --- a/vtk/src/layers/EGlyphLayer.cpp +++ b/vtk/src/layers/EGlyphLayer.cpp @@ -11,12 +11,36 @@ #include #include #include +#include #include #include "../CartographicTransformation.h" -#include "../advection/readdata.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() { diff --git a/vtk/src/layers/LGlyphLayer.cpp b/vtk/src/layers/LGlyphLayer.cpp index 9b0840a..845d628 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(std::unique_ptr advectionKernel) { +LGlyphLayer::LGlyphLayer() { this->ren = vtkSmartPointer::New(); this->ren->SetLayer(2); @@ -39,8 +39,6 @@ LGlyphLayer::LGlyphLayer(std::unique_ptr advectionKernel) { this->data = vtkSmartPointer::New(); this->data->SetPoints(this->points); - advector = std::move(advectionKernel); - auto camera = createNormalisedCamera(); ren->SetActiveCamera(camera); @@ -90,7 +88,7 @@ void LGlyphLayer::updateData(int t) { double point[3]; for (vtkIdType n = 0; n < this->points->GetNumberOfPoints(); n++) { this->points->GetPoint(n, point); - auto [yNew, xNew] = advector->advect(t, point[1], point[0]); + auto [xNew, yNew] = advect(n, point[0], point[1]); 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 1a8dc0f..22993b4 100644 --- a/vtk/src/layers/LGlyphLayer.h +++ b/vtk/src/layers/LGlyphLayer.h @@ -2,7 +2,6 @@ #define LGLYPHLAYER_H #include "Layer.h" -#include "../advection/AdvectionKernel.h" #include "../commands/SpawnPointCallback.h" #include #include @@ -14,14 +13,13 @@ class LGlyphLayer : public Layer { private: vtkSmartPointer points; vtkSmartPointer data; - std::unique_ptr advector; public: /** Constructor. */ - LGlyphLayer(std::unique_ptr advectionKernel); + LGlyphLayer(); /** 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 88c9bfa..fb89ea4 100644 --- a/vtk/src/main.cpp +++ b/vtk/src/main.cpp @@ -7,24 +7,18 @@ #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() { - shared_ptr uvGrid = std::make_shared(); - auto kernelRK4 = make_unique(uvGrid); - - auto l = new LGlyphLayer(move(kernelRK4)); -// l->spoofPoints(); + auto l = new LGlyphLayer(); + l->spoofPoints(); Program *program = new Program(); program->addLayer(new BackgroundImage("../../../../data/map_661-661.png"));