Added samples

CMakeLists.txt

@@ -3,15 +3,20 @@ cmake_minimum_required(VERSION 2.6)
 project(CosmoToolbox)
 
 find_path(NETCDF_INCLUDE_PATH NAMES netcdf.h)
+find_path(GSL_INCLUDE_PATH NAMES gsl/gsl_blas.h)
 
 find_library(NETCDF_LIBRARY netcdf)
 find_library(NETCDFCPP_LIBRARY netcdf_c++)
+find_library(GSL_LIBRARY gsl)
+find_library(GSLCBLAS_LIBRARY gslcblas)
 
 
 
 include(FindPackageHandleStandardArgs)
 set(NETCDF_FIND_REQUIRED TRUE)
+set(GSL_FIND_REQUIRED TRUE)
 FIND_PACKAGE_HANDLE_STANDARD_ARGS(NetCDF DEFAULT_MSG NETCDF_LIBRARY NETCDFCPP_LIBRARY NETCDF_INCLUDE_PATH)
+FIND_PACKAGE_HANDLE_STANDARD_ARGS(GSL DEFAULT_MSG GSL_LIBRARY GSLCBLAS_LIBRARY GSL_INCLUDE_PATH)
 
 
 # CPACK Configuration
@@ -27,5 +32,6 @@ SET(CPACK_SOURCE_IGNORE_FILES
 "/CVS/;/\\\\.git/;/\\\\.svn/;\\\\.swp$;\\\\.#;/#;.*~;cscope.*;/CMakeFiles/;.*\\\\.cmake;Makefile")
 
 add_subdirectory(src)
+add_subdirectory(sample)
 
 include(CPack)

sample/CMakeLists.txt

@@ -0,0 +1,11 @@
+SET(tolink ${CMAKE_BINARY_DIR}/src/libCosmoTool.so ${GSL_LIBRARY} ${GSLCBLAS_LIBRARY})
+include_directories(${CMAKE_SOURCE_DIR}/src)
+
+add_executable(testBQueue testBQueue.cpp)
+target_link_libraries(testBQueue ${tolink})
+
+add_executable(testInterpolate testInterpolate.cpp)
+target_link_libraries(testInterpolate ${tolink})
+
+add_executable(testSmooth testSmooth.cpp)
+target_link_libraries(testSmooth ${tolink})

sample/testBQueue.cpp

@@ -5,7 +5,7 @@ using namespace std;
 
 int main(int argc, char **argv)
 {
-  CosmoTool::BoundedQueue<int,int,100000> bq(4);
+  CosmoTool::BoundedQueue<int,int> bq(4, 100000.);
   
   for (int i = 10; i >= 0; i--)
     {

sample/testDelaunay.cpp

@@ -0,0 +1,28 @@
+#include <iostream>
+#include "dinterpolate.hpp"
+
+#define NX 10
+#define NY 10
+
+using namespace std;
+using namespace CosmoTool;
+
+typedef DelaunayInterpolate<double,double,2> myTriangle;
+
+int main()
+{
+  myTriangle::CoordType pos[] = { {0,0}, {1,0}, {0,1}, {1, 1} } ;
+  double vals[] = { 0, 1, 1, 0 };
+  uint32_t simplex[] = { 0, 1, 2, 3, 1, 2 };
+
+  myTriangle t(&pos[0], &vals[0], &simplex[0], 4, 2);
+
+  for (uint32_t iy = 0; iy <= NY; iy++) {
+    for (uint32_t ix = 0; ix <= NX; ix++) {
+      myTriangle::CoordType inter = { ix *1.0/ NX, iy *1.0/NY };
+      cout << inter[1] << " " << inter[0] << " " << t.computeValue(inter) << endl;
+    }
+    cout << endl;
+  }
+  return 0;
+}

sample/testDelaunay2.cpp

@@ -0,0 +1,26 @@
+#include <iostream>
+#include "dinterpolate.hpp"
+
+#define NX 10
+#define NY 10
+
+using namespace std;
+using namespace CosmoTool;
+
+typedef DelaunayInterpolate<double,double,2> myTriangle;
+
+int main()
+{
+
+
+  myTriangle t(&pos[0], &vals[0], &simplex[0], 4, 2);
+
+  for (uint32_t iy = 0; iy <= NY; iy++) {
+    for (uint32_t ix = 0; ix <= NX; ix++) {
+      myTriangle::CoordType inter = { ix *1.0/ NX, iy *1.0/NY };
+      cout << inter[1] << " " << inter[0] << " " << t.computeValue(inter) << endl;
+    }
+    cout << endl;
+  }
+  return 0;
+}

sample/testInterpolate.cpp

@@ -0,0 +1,30 @@
+#include <iostream>
+#include "interpolate3d.hpp"
+
+using namespace std;
+using namespace CosmoTool;
+
+int main()
+{
+   VectorField<float,2> *vectors = new VectorField<float,2>[8];
+
+   for (int i = 0; i < 2; i++)
+     for (int j = 0; j < 2; j++)
+       for (int k = 0; k < 2; k++)
+         {
+            int idx = i + 2*(j + 2*k);
+            vectors[idx].vec[0] = i;
+            vectors[idx].vec[1] = j;
+         }
+
+   GridSampler<VectorField<float,2> > sampler(vectors, 2, 2, 2, 2);
+   Interpolate3D<GridSampler<VectorField<float,2> > > inter(sampler);
+
+   VectorField<float,2> v = inter.get(0.5,0.5,0.5);
+   VectorField<float,2> v2 = inter.get(1.5,1.5,1.5);
+
+   cout << v.vec[0] << " " << v.vec[1] << endl;
+   cout << v2.vec[0] << " " << v2.vec[1] << endl;
+
+   return 0;
+}

sample/testSmooth.cpp

@@ -0,0 +1,68 @@
+#include <iostream>
+#include "sphSmooth.hpp"
+#include "yorick.hpp"
+#include "mykdtree.hpp"
+
+using namespace std;
+using namespace CosmoTool;
+
+#define NX 1024
+#define ND 2
+
+typedef SPHSmooth<char,ND> MySmooth;
+typedef MySmooth::SPHTree MyTree;
+typedef MySmooth::SPHNode MyNode;
+typedef MySmooth::SPHCell MyCell;
+
+double unit_fun(const char& c)
+{
+  return 1.0;
+}
+
+int main()
+{
+  uint32_t Ncells = 10000;
+  MyCell *cells = new MyCell[Ncells];
+  
+  for (int i = 0; i < Ncells; i++)
+    {
+      cells[i].active = true;
+      for (int l = 0; l < ND; l++)
+	cells[i].coord[l] = drand48();
+      cells[i].val.weight = 0;
+    }
+  
+  MyTree tree(cells, Ncells);
+  MySmooth smooth(&tree, 16);
+
+  for (uint32_t iy = 0; iy < NX; iy++)
+    {
+      cout << "iy=" << iy << endl;
+      for (uint32_t ix = 0; ix < NX; ix++)
+        {
+           MyTree::coords c = { 1.0*ix/NX, 1.0*iy/NX };
+           smooth.fetchNeighbours(c);
+           smooth.addGridSite(c);
+        }
+    }
+
+  
+  uint32_t dims[] = { NX, NX };
+  ProgressiveOutput<ComputePrecision> out = 
+    ProgressiveOutput<ComputePrecision>::saveArrayProgressive("out.nc", dims, 2);
+  for (uint32_t iy = 0; iy < NX; iy++)
+    {
+      cout << "iy=" << iy << endl;
+      for (uint32_t ix = 0; ix < NX; ix++)
+        {
+           MyTree::coords c = { 1.0*ix/NX, 1.0*iy/NX };
+           smooth.fetchNeighbours(c);
+           out.put(smooth.computeSmoothedValue(c, unit_fun));
+
+	  
+        }
+    }
+
+
+  return 0;
+}