py_nmie.cc 5.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122
  1. //**********************************************************************************//
  2. // Copyright (C) 2009-2017 Ovidio Pena <ovidio@bytesfall.com> //
  3. // Copyright (C) 2013-2017 Konstantin Ladutenko <kostyfisik@gmail.com> //
  4. // //
  5. // This file is part of scattnlay //
  6. // //
  7. // This program is free software: you can redistribute it and/or modify //
  8. // it under the terms of the GNU General Public License as published by //
  9. // the Free Software Foundation, either version 3 of the License, or //
  10. // (at your option) any later version. //
  11. // //
  12. // This program is distributed in the hope that it will be useful, //
  13. // but WITHOUT ANY WARRANTY; without even the implied warranty of //
  14. // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //
  15. // GNU General Public License for more details. //
  16. // //
  17. // The only additional remark is that we expect that all publications //
  18. // describing work using this software, or all commercial products //
  19. // using it, cite at least one of the following references: //
  20. // [1] O. Pena and U. Pal, "Scattering of electromagnetic radiation by //
  21. // a multilayered sphere," Computer Physics Communications, //
  22. // vol. 180, Nov. 2009, pp. 2348-2354. //
  23. // [2] K. Ladutenko, U. Pal, A. Rivera, and O. Pena-Rodriguez, "Mie //
  24. // calculation of electromagnetic near-field for a multilayered //
  25. // sphere," Computer Physics Communications, vol. 214, May 2017, //
  26. // pp. 225-230. //
  27. // //
  28. // You should have received a copy of the GNU General Public License //
  29. // along with this program. If not, see <http://www.gnu.org/licenses/>. //
  30. //**********************************************************************************//
  31. #include <math.h>
  32. #include <stdlib.h>
  33. #include <stdio.h>
  34. #include "nmie.hpp"
  35. #include "py_nmie.h"
  36. // Same as ScattCoeffs in 'nmie.h' but uses double arrays to return the results (useful for python).
  37. // This is a workaround because I have not been able to return the results using
  38. // std::vector<std::complex<double> >
  39. int ScattCoeffs(const unsigned int L, const int pl, std::vector<double>& x, std::vector<std::complex<double> >& m,
  40. const int nmax, double anr[], double ani[], double bnr[], double bni[]) {
  41. int i, result;
  42. std::vector<std::complex<double> > an, bn;
  43. an.resize(nmax);
  44. bn.resize(nmax);
  45. result = nmie::ScattCoeffs(L, pl, x, m, nmax, an, bn);
  46. for (i = 0; i < result; i++) {
  47. anr[i] = an[i].real();
  48. ani[i] = an[i].imag();
  49. bnr[i] = bn[i].real();
  50. bni[i] = bn[i].imag();
  51. }
  52. return result;
  53. }
  54. // Same as nMie in 'nmie.h' but uses double arrays to return the results (useful for python).
  55. // This is a workaround because I have not been able to return the results using
  56. // std::vector<std::complex<double> >
  57. int nMie(const int L, const int pl, std::vector<double>& x, std::vector<std::complex<double> >& m,
  58. const int nTheta, std::vector<double>& Theta, const int nmax,
  59. double *Qext, double *Qsca, double *Qabs, double *Qbk, double *Qpr, double *g, double *Albedo,
  60. double S1r[], double S1i[], double S2r[], double S2i[]) {
  61. int i, result;
  62. std::vector<std::complex<double> > S1, S2;
  63. S1.resize(nTheta);
  64. S2.resize(nTheta);
  65. result = nmie::nMie(L, pl, x, m, nTheta, Theta, nmax, Qext, Qsca, Qabs, Qbk, Qpr, g, Albedo, S1, S2);
  66. for (i = 0; i < nTheta; i++) {
  67. S1r[i] = S1[i].real();
  68. S1i[i] = S1[i].imag();
  69. S2r[i] = S2[i].real();
  70. S2i[i] = S2[i].imag();
  71. }
  72. return result;
  73. }
  74. // Same as nField in 'nmie.h' but uses double arrays to return the results (useful for python).
  75. // This is a workaround because I have not been able to return the results using
  76. // std::vector<std::complex<double> >
  77. int nField(const int L, const int pl, std::vector<double>& x, std::vector<std::complex<double> >& m, const int nmax,
  78. const int nCoords, std::vector<double>& Xp, std::vector<double>& Yp, std::vector<double>& Zp,
  79. double Erx[], double Ery[], double Erz[], double Eix[], double Eiy[], double Eiz[],
  80. double Hrx[], double Hry[], double Hrz[], double Hix[], double Hiy[], double Hiz[]) {
  81. int i, result;
  82. std::vector<std::vector<std::complex<double> > > E, H;
  83. E.resize(nCoords);
  84. H.resize(nCoords);
  85. for (i = 0; i < nCoords; i++) {
  86. E[i].resize(3);
  87. H[i].resize(3);
  88. }
  89. result = nmie::nField(L, pl, x, m, nmax, nCoords, Xp, Yp, Zp, E, H);
  90. for (i = 0; i < nCoords; i++) {
  91. Erx[i] = E[i][0].real();
  92. Ery[i] = E[i][1].real();
  93. Erz[i] = E[i][2].real();
  94. Eix[i] = E[i][0].imag();
  95. Eiy[i] = E[i][1].imag();
  96. Eiz[i] = E[i][2].imag();
  97. Hrx[i] = H[i][0].real();
  98. Hry[i] = H[i][1].real();
  99. Hrz[i] = H[i][2].real();
  100. Hix[i] = H[i][0].imag();
  101. Hiy[i] = H[i][1].imag();
  102. Hiz[i] = H[i][2].imag();
  103. }
  104. return result;
  105. }