Browse Source

scattcoeffs works for single set

Konstantin Ladutenko 6 years ago
parent
commit
662b0b1b3c
2 changed files with 31 additions and 28 deletions
  1. 6 1
      examples/scattcoeffs.py
  2. 25 27
      src/nmie.cc

+ 6 - 1
examples/scattcoeffs.py

@@ -51,6 +51,8 @@
 from scattnlay import scattcoeffs
 import numpy as np
 
+import example
+
 size = np.arange(0.25, 100.25, 0.25)
 
 x = np.ones((len(size), 5), dtype = np.float64)
@@ -68,6 +70,9 @@ m[:, 3] *= 2.8 + 0.2j
 m[:, 4] *= 1.5 + 0.4j
 
 terms, an, bn = scattcoeffs(x, m, 105)
+terms1, an1, bn1 = example.scattcoeffs(x[0,:], m[0,:])
+print(terms)
+print(terms1)
 
 result = np.vstack((x[:, 4], an[:, 0].real, an[:, 0].imag, an[:, 1].real, an[:, 1].imag, an[:, 2].real, an[:, 2].imag,
                              bn[:, 0].real, bn[:, 0].imag, bn[:, 1].real, bn[:, 1].imag, bn[:, 2].real, bn[:, 2].imag)).transpose()
@@ -91,5 +96,5 @@ try:
     plt.show()
 finally:
     np.savetxt("scattcoeffs.txt", result, fmt = "%.5f")
-    print result
+    print( result)
 

+ 25 - 27
src/nmie.cc

@@ -431,42 +431,40 @@ namespace nmie {
   // int ScattCoeffs(const unsigned int L, const int pl, std::vector<double>& x, std::vector<std::complex<double> >& m,
   //                 const int nmax, std::vector<std::complex<double> >& an, std::vector<std::complex<double> >& bn) {
 
-    
-std::complex<double>
-test_stdvec(const  std::complex<double> a, const std::vector< std::complex<double> > input, std::vector< std::complex<double> >& output) {
-  for (unsigned int i = 0; i < input.size();++i){    
-    output[i] = input[i]*a;
-  }
-  return a;
-}
 
-std::tuple<
-  py::array_t< std::complex<double>, py::array::c_style | py::array::forcecast>,
- std::complex<double>
-  >
-py_test_stdvec(const  std::complex<double> a,
-                   py::array_t< std::complex<double>, py::array::c_style | py::array::forcecast> py_input) {
-  std::vector< std::complex<double> > c_input(py_input.size());
-  std::vector< std::complex<double> > c_output(py_input.size());
+py::tuple py_ScattCoeffs(
+               py::array_t<double, py::array::c_style | py::array::forcecast> py_x,
+               py::array_t< std::complex<double>, py::array::c_style | py::array::forcecast> py_m// ,
+               // const int nmax=-1, const int pl=-1
+               ) {
   
-  std::memcpy(c_input.data(), py_input.data(), py_input.size()*sizeof( std::complex<double>));
-
-   std::complex<double> c_a = 0;
-  c_a = test_stdvec(a, c_input, c_output);
-    
+  std::vector<double> c_x(py_x.size());
+  std::vector< std::complex<double> > c_m(py_m.size());
+  int L = py_x.size();
+  int nmax = -1, pl = -1;
+  
+  std::memcpy(c_x.data(), py_x.data(), py_x.size()*sizeof(double));
+  std::memcpy(c_m.data(), py_m.data(), py_m.size()*sizeof( std::complex<double>));
 
-  auto result        = py::array_t< std::complex<double>>(c_output.size());
-  auto result_buffer = result.request();
-   std::complex<double> *result_ptr    = ( std::complex<double> *) result_buffer.ptr;
+  int terms = 0;
+  std::vector<std::complex<double> > c_an, c_bn;
+  terms = nmie::ScattCoeffs( L, pl, c_x, c_m, nmax, c_an, c_bn);
 
-  std::memcpy(result_ptr,c_output.data(),c_output.size()*sizeof( std::complex<double>));
+  auto py_an = py::array_t< std::complex<double>>(c_an.size());
+  auto py_bn = py::array_t< std::complex<double>>(c_bn.size());
+  auto py_an_buffer = py_an.request();
+  auto py_bn_buffer = py_bn.request();
+  std::complex<double> *py_an_ptr    = ( std::complex<double> *) py_an_buffer.ptr;
+  std::complex<double> *py_bn_ptr    = ( std::complex<double> *) py_bn_buffer.ptr;
 
-  return std::make_tuple(result, c_a);
+  std::memcpy(py_an_ptr,c_an.data(),c_an.size()*sizeof( std::complex<double>));
+  std::memcpy(py_bn_ptr,c_bn.data(),c_bn.size()*sizeof( std::complex<double>));
 
+  return py::make_tuple(terms, py_an, py_bn);
 }
 
 PYBIND11_MODULE(example, m) {
     m.doc() = "pybind11 example plugin"; // optional module docstring
 
-    m.def("test_stdvec", &py_test_stdvec, "Multipy vec by a const");
+    m.def("scattcoeffs", &py_ScattCoeffs, "test");
 }