ARTS 2.5.0 (git: 9ee3ac6c)
tmatrix.cc
Go to the documentation of this file.
1/* Copyright (C) 2013 Oliver Lemke
2 *
3 * This program is free software: you can redistribute it and/or modify
4 * it under the terms of the GNU General Public License as published by
5 * the Free Software Foundation, either version 3 of the License, or
6 * (at your option) any later version.
7 *
8 * This program is distributed in the hope that it will be useful,
9 * but WITHOUT ANY WARRANTY; without even the implied warranty of
10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 * GNU General Public License for more details.
12 *
13 * You should have received a copy of the GNU General Public License
14 * along with this program. If not, see <http://www.gnu.org/licenses/>.
15 *
16 */
17
26#include "tmatrix.h"
27#include <cmath>
28#include <cstring>
29#include <stdexcept>
30#include "matpack_complex.h"
31#include "math_funcs.h"
32#include "matpackI.h"
33#include "messages.h"
34#include "optproperties.h"
35
36void calc_phamat(Matrix& z,
37 const Index& nmax,
38 const Numeric& lam,
39 const Numeric& thet0,
40 const Numeric& thet,
41 const Numeric& phi0,
42 const Numeric& phi,
43 const Numeric& beta,
44 const Numeric& alpha);
45
47 const Complex& s11,
48 const Complex& s12,
49 const Complex& s21,
50 const Complex& s22);
51
53 const Index& nmax,
54 const Numeric& lam,
55 const Numeric& thet0,
56 const Numeric& thet,
57 const Numeric& phi0,
58 const Numeric& phi,
59 const Numeric& beta,
60 const Numeric& alpha_1,
61 const Numeric& alpha_2);
62
64 const Index& nmax,
65 const Numeric& lam,
66 const Numeric& thet0,
67 const Numeric& thet,
68 const Numeric& phi0,
69 const Numeric& phi,
70 const Numeric& beta,
71 const Numeric& alpha_1,
72 const Numeric& alpha_2);
73
75 const Index& nmax,
76 const Numeric& lam,
77 const Numeric& thet0_1,
78 const Numeric& thet0_2,
79 const Numeric& thet,
80 const Numeric& phi0,
81 const Numeric& phi_1,
82 const Numeric& phi_2,
83 const Numeric& beta,
84 const Numeric& alpha);
85
87 const Index& nmax,
88 const Numeric& lam,
89 const Numeric& thet0_1,
90 const Numeric& thet0_2,
91 const Numeric& thet,
92 const Numeric& phi0,
93 const Numeric& phi_1,
94 const Numeric& phi_2,
95 const Numeric& beta,
96 const Numeric& alpha_1,
97 const Numeric& alpha_2);
98
99#ifdef ENABLE_TMATRIX
100extern "C" {
101#endif
209void tmd_(const Numeric& rat,
210 const Index& ndistr,
211 const Numeric& axmax,
212 const Index& npnax,
213 const Numeric& b,
214 const Numeric& gam,
215 const Index& nkmax,
216 const Numeric& eps,
217 const Index& np,
218 const Numeric& lam,
219 const Numeric& mrr,
220 const Numeric& mri,
221 const Numeric& ddelt,
222 const Index& npna,
223 const Index& ndgs,
224 const Numeric& r1rat,
225 const Numeric& r2rat,
226 const Index& quiet,
227 Numeric& reff, // OUT
228 Numeric& veff, // OUT
229 Numeric& cext, // OUT
230 Numeric& csca, // OUT
231 Numeric& walb, // OUT
232 Numeric& asymm, // OUT
233 Numeric* f11, // Array npna
234 Numeric* f22, // Array npna
235 Numeric* f33, // Array npna
236 Numeric* f44, // Array npna
237 Numeric* f12, // Array npna
238 Numeric* f34, // Array npna
239 char* errmsg);
240
274void tmatrix_(const Numeric& rat,
275 const Numeric& axi,
276 const Index& np,
277 const Numeric& lam,
278 const Numeric& eps,
279 const Numeric& mrr,
280 const Numeric& mri,
281 const Numeric& ddelt,
282 const Index& quiet,
283 Index& nmax,
284 Numeric& csca,
285 Numeric& cext,
286 char* errmsg);
287
311void ampl_(const Index& nmax,
312 const Numeric& lam,
313 const Numeric& thet0,
314 const Numeric& thet,
315 const Numeric& phi0,
316 const Numeric& phi,
317 const Numeric& alpha,
318 const Numeric& beta,
319 Complex& s11,
320 Complex& s12,
321 Complex& s21,
322 Complex& s22);
323
332void avgtmatrix_(const Index& nmax);
333#ifdef ENABLE_TMATRIX
334}
335#endif
336
337// Define dummy functions that throw runtime errors if ARTS is
338// compiled without T-Matrix support.
339#ifndef ENABLE_TMATRIX
340
341// T-matrix code for randomly oriented nonspherical particles.
342void tmd_(const Numeric&,
343 const Index&,
344 const Numeric&,
345 const Index&,
346 const Numeric&,
347 const Numeric&,
348 const Index&,
349 const Numeric&,
350 const Index&,
351 const Numeric&,
352 const Numeric&,
353 const Numeric&,
354 const Numeric&,
355 const Index&,
356 const Index&,
357 const Numeric&,
358 const Numeric&,
359 const Index&,
360 Numeric&,
361 Numeric&,
362 Numeric&,
363 Numeric&,
364 Numeric&,
365 Numeric&,
366 Numeric*,
367 Numeric*,
368 Numeric*,
369 Numeric*,
370 Numeric*,
371 Numeric*,
372 char*) {
373 throw std::runtime_error(
374 "This version of ARTS was compiled without T-Matrix support.");
375}
376
377// T-matrix code for nonspherical particles in a fixed orientation
378void tmatrix_(const Numeric&,
379 const Numeric&,
380 const Index&,
381 const Numeric&,
382 const Numeric&,
383 const Numeric&,
384 const Numeric&,
385 const Numeric&,
386 const Index&,
387 Index&,
388 Numeric&,
389 Numeric&,
390 char*) {
391 throw std::runtime_error(
392 "This version of ARTS was compiled without T-Matrix support.");
393}
394
395// T-matrix code for nonspherical particles in a fixed orientation
396void ampl_(const Index&,
397 const Numeric&,
398 const Numeric&,
399 const Numeric&,
400 const Numeric&,
401 const Numeric&,
402 const Numeric&,
403 const Numeric&,
404 Complex&,
405 Complex&,
406 Complex&,
407 Complex&) {
408 throw std::runtime_error(
409 "This version of ARTS was compiled without T-Matrix support.");
410}
411
412void avgtmatrix_(const Index&) {
413 throw std::runtime_error(
414 "This version of ARTS was compiled without T-Matrix support.");
415}
416
417#endif
418
420 const Index& nmax,
421 const Numeric& lam,
422 const Numeric& thet0,
423 const Numeric& thet,
424 const Numeric& phi0,
425 const Numeric& phi,
426 const Numeric& beta,
427 const Numeric& alpha) {
428 Complex s11;
429 Complex s12;
430 Complex s21;
431 Complex s22;
432 ampl_(nmax, lam, thet0, thet, phi0, phi, alpha, beta, s11, s12, s21, s22);
433
434 ampmat_to_phamat(z, s11, s12, s21, s22);
435}
436
451 const Complex& s11,
452 const Complex& s12,
453 const Complex& s21,
454 const Complex& s22) {
455 z.resize(4, 4);
456 z(0, 0) = 0.5 * (s11 * conj(s11) + s12 * conj(s12) + s21 * conj(s21) +
457 s22 * conj(s22))
458 .real();
459 z(0, 1) = 0.5 * (s11 * conj(s11) - s12 * conj(s12) + s21 * conj(s21) -
460 s22 * conj(s22))
461 .real();
462 z(0, 2) = (-s11 * conj(s12) - s22 * conj(s21)).real();
463 z(0, 3) = (Complex(0., 1.) * (s11 * conj(s12) - s22 * conj(s21))).real();
464
465 z(1, 0) = 0.5 * (s11 * conj(s11) + s12 * conj(s12) - s21 * conj(s21) -
466 s22 * conj(s22))
467 .real();
468 z(1, 1) = 0.5 * (s11 * conj(s11) - s12 * conj(s12) - s21 * conj(s21) +
469 s22 * conj(s22))
470 .real();
471 z(1, 2) = (-s11 * conj(s12) + s22 * conj(s21)).real();
472 z(1, 3) = (Complex(0., 1.) * (s11 * conj(s12) + s22 * conj(s21))).real();
473
474 z(2, 0) = (-s11 * conj(s21) - s22 * conj(s12)).real();
475 z(2, 1) = (-s11 * conj(s21) + s22 * conj(s12)).real();
476 z(2, 2) = (s11 * conj(s22) + s12 * conj(s21)).real();
477 z(2, 3) = (Complex(0., -1.) * (s11 * conj(s22) + s21 * conj(s12))).real();
478
479 z(3, 0) = (Complex(0., 1.) * (s21 * conj(s11) + s22 * conj(s12))).real();
480 z(3, 1) = (Complex(0., 1.) * (s21 * conj(s11) - s22 * conj(s12))).real();
481 z(3, 2) = (Complex(0., -1.) * (s22 * conj(s11) - s12 * conj(s21))).real();
482 z(3, 3) = (s22 * conj(s11) - s12 * conj(s21)).real();
483}
484
485static const Numeric GaussLeg6[][3] = {{0.23861918, 0.66120939, 0.93246951},
486 {0.46791393, 0.36076157, 0.17132449}};
487
488static const Numeric GaussLeg10[][5] = {
489 {0.14887434, 0.43339539, 0.67940957, 0.86506337, 0.97390653},
490 {0.29552422, 0.26926672, 0.21908636, 0.14945135, 0.06667134}};
491
511 const Index& nmax,
512 const Numeric& lam,
513 const Numeric& thet0,
514 const Numeric& thet,
515 const Numeric& phi0,
516 const Numeric& phi,
517 const Numeric& beta,
518 const Numeric& alpha_1,
519 const Numeric& alpha_2) {
520 const Numeric c = 0.5 * (alpha_2 + alpha_1);
521 const Numeric m = 0.5 * (alpha_2 - alpha_1);
522
523 phamat.resize(4, 4);
524 phamat = 0.;
525 Matrix z;
526
527 for (Index i = 0; i < 5; ++i) {
528 const Numeric abscissa = GaussLeg10[0][i];
529 const Numeric weight = GaussLeg10[1][i];
530
531 calc_phamat(z, nmax, lam, thet0, thet, phi0, phi, beta, c + m * abscissa);
532 z *= weight;
533 phamat += z;
534
535 calc_phamat(z, nmax, lam, thet0, thet, phi0, phi, beta, c - m * abscissa);
536 z *= weight;
537 phamat += z;
538 }
539 phamat *= m;
540}
541
561 const Index& nmax,
562 const Numeric& lam,
563 const Numeric& thet0,
564 const Numeric& thet,
565 const Numeric& phi0,
566 const Numeric& phi,
567 const Numeric& beta,
568 const Numeric& alpha_1,
569 const Numeric& alpha_2) {
570 const Numeric c = 0.5 * (alpha_2 + alpha_1);
571 const Numeric m = 0.5 * (alpha_2 - alpha_1);
572
573 phamat.resize(4, 4);
574 phamat = 0.;
575 Matrix z;
576
577 for (Index i = 0; i < 3; ++i) {
578 const Numeric abscissa = GaussLeg6[0][i];
579 const Numeric weight = GaussLeg6[1][i];
580
581 calc_phamat(z, nmax, lam, thet0, thet, phi0, phi, beta, c + m * abscissa);
582 z *= weight;
583 phamat += z;
584
585 calc_phamat(z, nmax, lam, thet0, thet, phi0, phi, beta, c - m * abscissa);
586 z *= weight;
587 phamat += z;
588 }
589 phamat *= m;
590}
591
613 const Index& nmax,
614 const Numeric& lam,
615 const Numeric& thet0_1,
616 const Numeric& thet0_2,
617 const Numeric& thet,
618 const Numeric& phi0,
619 const Numeric& phi_1,
620 const Numeric& phi_2,
621 const Numeric& beta,
622 const Numeric& alpha) {
623 extern const Numeric PI;
624
625 phamat.resize(4, 4);
626 phamat = 0.;
627 Matrix z;
628
629 const Numeric c_thet0 = 0.5 * (thet0_2 + thet0_1);
630 const Numeric m_thet0 = 0.5 * (thet0_2 - thet0_1);
631 const Numeric c_phi = 0.5 * (phi_2 + phi_1);
632 const Numeric m_phi = 0.5 * (phi_2 - phi_1);
633
634 for (Index t = 0; t < 5; ++t) {
635 const Numeric abscissa_thet0 = GaussLeg10[0][t];
636 const Numeric weight_thet0 = GaussLeg10[1][t];
637
638 Matrix phamat_phi(4, 4, 0.);
639
640 for (Index p = 0; p < 5; ++p) {
641 const Numeric abscissa_phi = GaussLeg10[0][p];
642 const Numeric weight_phi = GaussLeg10[1][p];
643
644 Numeric this_thet0 = c_thet0 + m_thet0 * abscissa_thet0;
645 calc_phamat(z,
646 nmax,
647 lam,
648 this_thet0,
649 thet,
650 phi0,
651 c_phi + m_phi * abscissa_phi,
652 beta,
653 alpha);
654 z *= weight_phi * sin(this_thet0 * PI / 180.);
655 phamat_phi += z;
656
657 calc_phamat(z,
658 nmax,
659 lam,
660 this_thet0,
661 thet,
662 phi0,
663 c_phi - m_phi * abscissa_phi,
664 beta,
665 alpha);
666 z *= weight_phi * sin(this_thet0 * PI / 180.);
667 phamat_phi += z;
668
669 this_thet0 = c_thet0 - m_thet0 * abscissa_thet0;
670 calc_phamat(z,
671 nmax,
672 lam,
673 this_thet0,
674 thet,
675 phi0,
676 c_phi + m_phi * abscissa_phi,
677 beta,
678 alpha);
679 z *= weight_phi * sin(this_thet0 * PI / 180.);
680 phamat_phi += z;
681
682 calc_phamat(z,
683 nmax,
684 lam,
685 this_thet0,
686 thet,
687 phi0,
688 c_phi - m_phi * abscissa_phi,
689 beta,
690 alpha);
691 z *= weight_phi * sin(this_thet0 * PI / 180.);
692 phamat_phi += z;
693 }
694 phamat_phi *= m_phi * weight_thet0;
695 phamat += phamat_phi;
696 }
697
698 phamat *= m_thet0;
699}
700
724 const Index& nmax,
725 const Numeric& lam,
726 const Numeric& thet0_1,
727 const Numeric& thet0_2,
728 const Numeric& thet,
729 const Numeric& phi0,
730 const Numeric& phi_1,
731 const Numeric& phi_2,
732 const Numeric& beta,
733 const Numeric& alpha_1,
734 const Numeric& alpha_2) {
735 extern const Numeric PI;
736
737 Matrix z;
738 Matrix phamat_phi(4, 4);
739
740 const Numeric c_thet0 = 0.5 * (thet0_2 + thet0_1);
741 const Numeric m_thet0 = 0.5 * (thet0_2 - thet0_1);
742 const Numeric c_phi = 0.5 * (phi_2 + phi_1);
743 const Numeric m_phi = 0.5 * (phi_2 - phi_1);
744
745 phamat.resize(4, 4);
746 phamat = 0.;
747
748 for (Index t = 0; t < 3; ++t) {
749 const Numeric abscissa_thet0 = GaussLeg6[0][t];
750 const Numeric weight_thet0 = GaussLeg6[1][t];
751
752 phamat_phi = 0.;
753
754 for (Index p = 0; p < 3; ++p) {
755 const Numeric abscissa_phi = GaussLeg6[0][p];
756 const Numeric weight_phi = GaussLeg6[1][p];
757
758 Numeric this_thet0 = c_thet0 + m_thet0 * abscissa_thet0;
760 nmax,
761 lam,
762 this_thet0,
763 thet,
764 phi0,
765 c_phi + m_phi * abscissa_phi,
766 beta,
767 alpha_1,
768 alpha_2);
769 z *= weight_phi * sin(this_thet0 * PI / 180.);
770 phamat_phi += z;
771
773 nmax,
774 lam,
775 this_thet0,
776 thet,
777 phi0,
778 c_phi - m_phi * abscissa_phi,
779 beta,
780 alpha_1,
781 alpha_2);
782 z *= weight_phi * sin(this_thet0 * PI / 180.);
783 phamat_phi += z;
784
785 this_thet0 = c_thet0 - m_thet0 * abscissa_thet0;
787 nmax,
788 lam,
789 this_thet0,
790 thet,
791 phi0,
792 c_phi + m_phi * abscissa_phi,
793 beta,
794 alpha_1,
795 alpha_2);
796 z *= weight_phi * sin(this_thet0 * PI / 180.);
797 phamat_phi += z;
798
800 nmax,
801 lam,
802 this_thet0,
803 thet,
804 phi0,
805 c_phi - m_phi * abscissa_phi,
806 beta,
807 alpha_1,
808 alpha_2);
809 z *= weight_phi * sin(this_thet0 * PI / 180.);
810 phamat_phi += z;
811 }
812 phamat_phi *= m_phi * weight_thet0;
813 phamat += phamat_phi;
814 }
815
816 phamat *= m_thet0;
817}
818
844 Numeric& csca,
845 Vector& f11,
846 Vector& f22,
847 Vector& f33,
848 Vector& f44,
849 Vector& f12,
850 Vector& f34,
851 const Numeric equiv_radius,
852 const Numeric aspect_ratio,
853 const Index np,
854 const Numeric lam,
855 const Numeric ref_index_real,
856 const Numeric ref_index_imag,
857 const Numeric precision,
858 const Index nza,
859 const Index ndgs,
860 const Index quiet = 1) {
861 Numeric reff;
862 Numeric veff;
863 Numeric walb;
864 Numeric asymm;
865
866 char errmsg[1024] = "";
867
868 f11.resize(nza);
869 f11 = NAN;
870 f22.resize(nza);
871 f22 = NAN;
872 f33.resize(nza);
873 f33 = NAN;
874 f44.resize(nza);
875 f44 = NAN;
876 f12.resize(nza);
877 f12 = NAN;
878 f34.resize(nza);
879 f34 = NAN;
880
881 // It is necessary to make sure that the Fortran code is not
882 // called from different threads at the same time. The common
883 // blocks are not threadsafe.
884#pragma omp critical(tmatrix_code)
885 tmd_(1.0,
886 4,
887 equiv_radius,
888 1,
889 0.1,
890 1.0,
891 -1,
892 aspect_ratio,
893 np,
894 lam,
895 ref_index_real,
896 ref_index_imag,
897 precision,
898 nza,
899 ndgs,
900 0.9999999,
901 1.0000001,
902 quiet,
903 reff,
904 veff,
905 cext,
906 csca,
907 walb,
908 asymm,
909 f11.get_c_array(),
910 f22.get_c_array(),
911 f33.get_c_array(),
912 f44.get_c_array(),
913 f12.get_c_array(),
914 f34.get_c_array(),
915 errmsg);
916
917 if (strlen(errmsg)) {
918 std::ostringstream os;
919 os << "T-Matrix code failed: " << errmsg;
920 throw std::runtime_error(os.str());
921 }
922}
923
942 Numeric& csca,
943 Index& nmax,
944 const Numeric equiv_radius,
945 const Numeric aspect_ratio,
946 const Index np,
947 const Numeric lam,
948 const Numeric ref_index_real,
949 const Numeric ref_index_imag,
950 const Numeric precision,
951 const Index quiet = 1) {
952 char errmsg[1024] = "";
953
954 // It is necessary to make sure that the Fortran code is not
955 // called from different threads at the same time. The common
956 // blocks are not threadsafe.
957#pragma omp critical(tmatrix_code)
958 tmatrix_(1.,
959 equiv_radius,
960 np,
961 lam,
962 aspect_ratio,
963 ref_index_real,
964 ref_index_imag,
965 precision,
966 quiet,
967 nmax,
968 csca,
969 cext,
970 errmsg);
971
972 if (strlen(errmsg)) {
973 std::ostringstream os;
974 os << "T-Matrix code failed: " << errmsg;
975 throw std::runtime_error(os.str());
976 }
977}
978
980 ConstMatrixView ref_index_real,
981 ConstMatrixView ref_index_imag,
982 const Numeric equiv_radius,
983 const Index np,
984 const Numeric aspect_ratio,
985 const Numeric precision,
986 const Index ndgs,
987 const Index robust,
988 const Index quiet) {
989 const Index nf = ssd.f_grid.nelem();
990 const Index nT = ssd.T_grid.nelem();
991
992 extern const Numeric PI;
993 extern const Numeric SPEED_OF_LIGHT;
994
995 if (ref_index_real.nrows() != nf)
996 throw std::runtime_error(
997 "Number of rows of refractive index real part must match ssd.f_grid.");
998 if (ref_index_real.ncols() != nT)
999 throw std::runtime_error(
1000 "Number of cols of refractive index real part must match ssd.T_grid.");
1001 if (ref_index_imag.nrows() != nf)
1002 throw std::runtime_error(
1003 "Number of rows of refractive index imaginary part must match ssd.f_grid.");
1004 if (ref_index_imag.ncols() != nT)
1005 throw std::runtime_error(
1006 "Number of cols of refractive index imaginary part must match ssd.T_grid.");
1007
1008 // The T-Matrix code needs wavelength
1009 Vector lam(nf, SPEED_OF_LIGHT);
1010 lam /= ssd.f_grid;
1011
1012 switch (ssd.ptype) {
1013 case PTYPE_TOTAL_RND: {
1014 const Index nza = ssd.za_grid.nelem();
1015
1016 ssd.pha_mat_data.resize(nf, nT, nza, 1, 1, 1, 6);
1017 ssd.ext_mat_data.resize(nf, nT, 1, 1, 1);
1018 ssd.abs_vec_data.resize(nf, nT, 1, 1, 1);
1019
1020 ssd.pha_mat_data = NAN;
1021 ssd.ext_mat_data = NAN;
1022 ssd.abs_vec_data = NAN;
1023
1024 // Output variables
1025 Numeric cext = NAN;
1026 Numeric csca = NAN;
1027 Vector f11;
1028 Vector f22;
1029 Vector f33;
1030 Vector f44;
1031 Vector f12;
1032 Vector f34;
1033 Matrix mono_pha_mat_data(nza, 6, NAN);
1034
1035 ostringstream os;
1036 os << "Calculation of SingleScatteringData properties failed for\n\n";
1037 bool anyfailed = false;
1038#pragma omp critical(tmatrix_ssp)
1039 for (Index f_index = 0; f_index < nf; ++f_index)
1040 for (Index T_index = 0; T_index < nT; ++T_index) {
1041 bool thisfailed = false;
1042 try {
1044 csca,
1045 f11,
1046 f22,
1047 f33,
1048 f44,
1049 f12,
1050 f34,
1051 equiv_radius,
1052 aspect_ratio,
1053 np,
1054 lam[f_index],
1055 ref_index_real(f_index, T_index),
1056 ref_index_imag(f_index, T_index),
1057 precision,
1058 nza,
1059 ndgs,
1060 quiet);
1061 } catch (const std::runtime_error& e) {
1062 //ostringstream os;
1063 //os << "Calculation of SingleScatteringData properties failed for\n"
1064 os << "f_grid[" << f_index << "] = " << ssd.f_grid[f_index] << "\n"
1065 << "T_grid[" << T_index << "] = " << ssd.T_grid[T_index]
1066 << "\n\n";
1067 //<< e.what();
1068 //throw std::runtime_error(os.str());
1069 thisfailed = true;
1070 anyfailed = true;
1071 cout << "\n\n";
1072 }
1073
1074 if (!thisfailed) {
1075 mono_pha_mat_data(joker, 0) = f11;
1076 mono_pha_mat_data(joker, 1) = f12;
1077 mono_pha_mat_data(joker, 2) = f22;
1078 mono_pha_mat_data(joker, 3) = f33;
1079 mono_pha_mat_data(joker, 4) = f34;
1080 mono_pha_mat_data(joker, 5) = f44;
1081
1082 mono_pha_mat_data *= csca / 4. / PI;
1083 ssd.pha_mat_data(f_index, T_index, joker, 0, 0, 0, joker) =
1084 mono_pha_mat_data;
1085
1086 ssd.ext_mat_data(f_index, T_index, 0, 0, 0) = cext;
1087 ssd.abs_vec_data(f_index, T_index, 0, 0, 0) = cext - csca;
1088 }
1089 }
1090 if (anyfailed)
1091 if (robust)
1092 cout << os.str();
1093 else
1094 throw std::runtime_error(os.str());
1095 else {
1096 os << "None\n";
1097 }
1098
1099 break;
1100 }
1101 case PTYPE_AZIMUTH_RND: {
1102 const Index nza = ssd.za_grid.nelem();
1103 const Index naa = ssd.aa_grid.nelem();
1104
1105 ssd.pha_mat_data.resize(nf, nT, nza, naa, nza, 1, 16);
1106 ssd.ext_mat_data.resize(nf, nT, nza, 1, 3);
1107 ssd.abs_vec_data.resize(nf, nT, nza, 1, 2);
1108
1109 ssd.ext_mat_data = NAN;
1110 ssd.pha_mat_data = NAN;
1111 ssd.abs_vec_data = NAN;
1112
1113 // Output variables
1114 Numeric cext = NAN;
1115 Numeric csca = NAN;
1116 Index nmax = -1;
1117
1118 Tensor5 csca_data(nf, nT, nza, 1, 2);
1119
1120#pragma omp critical(tmatrix_ssp)
1121 for (Index f_index = 0; f_index < nf; ++f_index) {
1122 const Numeric lam_f = lam[f_index];
1123
1124 for (Index T_index = 0; T_index < nT; ++T_index) {
1125 try {
1127 csca,
1128 nmax,
1129 equiv_radius,
1130 aspect_ratio,
1131 np,
1132 lam_f,
1133 ref_index_real(f_index, T_index),
1134 ref_index_imag(f_index, T_index),
1135 precision);
1136 } catch (const std::runtime_error& e) {
1137 ostringstream os;
1138 os << "Calculation of SingleScatteringData properties failed for\n"
1139 << "f_grid[" << f_index << "] = " << ssd.f_grid[f_index] << "\n"
1140 << "T_grid[" << T_index << "] = " << ssd.T_grid[T_index] << "\n"
1141 << e.what();
1142 throw std::runtime_error(os.str());
1143 }
1144
1145 Matrix phamat;
1146 for (Index za_scat_index = 0; za_scat_index < nza; ++za_scat_index)
1147 for (Index aa_index = 0; aa_index < naa; ++aa_index)
1148 for (Index za_inc_index = 0; za_inc_index < nza; ++za_inc_index) {
1149 if (aspect_ratio < 1.0) {
1150 // Phase matrix for prolate particles
1152 nmax,
1153 lam_f,
1154 ssd.za_grid[za_inc_index],
1155 ssd.za_grid[za_scat_index],
1156 0.0,
1157 ssd.aa_grid[aa_index],
1158 90.0,
1159 0.0,
1160 180.0);
1161 phamat /= 180.;
1162 } else {
1163 // Phase matrix for oblate particles
1164 calc_phamat(phamat,
1165 nmax,
1166 lam_f,
1167 ssd.za_grid[za_inc_index],
1168 ssd.za_grid[za_scat_index],
1169 0.0,
1170 ssd.aa_grid[aa_index],
1171 0.0,
1172 0.0);
1173 }
1174
1175 ssd.pha_mat_data(f_index,
1176 T_index,
1177 za_scat_index,
1178 aa_index,
1179 za_inc_index,
1180 0,
1181 Range(0, 4)) = phamat(0, joker);
1182 ssd.pha_mat_data(f_index,
1183 T_index,
1184 za_scat_index,
1185 aa_index,
1186 za_inc_index,
1187 0,
1188 Range(4, 4)) = phamat(1, joker);
1189 ssd.pha_mat_data(f_index,
1190 T_index,
1191 za_scat_index,
1192 aa_index,
1193 za_inc_index,
1194 0,
1195 Range(8, 4)) = phamat(2, joker);
1196 ssd.pha_mat_data(f_index,
1197 T_index,
1198 za_scat_index,
1199 aa_index,
1200 za_inc_index,
1201 0,
1202 Range(12, 4)) = phamat(3, joker);
1203 }
1204
1205 // Csca integral
1206 for (Index za_scat_index = 0; za_scat_index < nza; ++za_scat_index) {
1207 Matrix csca_integral;
1208 if (aspect_ratio < 1.0) {
1209 // Csca for prolate particles
1211 nmax,
1212 lam_f,
1213 0,
1214 180,
1215 ssd.za_grid[za_scat_index],
1216 0.,
1217 0.,
1218 180.,
1219 90.,
1220 0.,
1221 180.);
1222 csca_integral /= 180.;
1223 } else {
1224 // Csca for oblate particles
1225 integrate_phamat_theta0_phi10(csca_integral,
1226 nmax,
1227 lam_f,
1228 0,
1229 180,
1230 ssd.za_grid[za_scat_index],
1231 0.,
1232 0.,
1233 180,
1234 0.,
1235 0.);
1236 }
1237 csca_data(f_index, T_index, za_scat_index, 0, joker) =
1238 csca_integral(Range(0, 2), 0);
1239 }
1240
1241 // Extinction matrix
1242 if (aspect_ratio < 1.0) {
1243 // Average T-Matrix for prolate particles
1244 avgtmatrix_(nmax);
1245 }
1246
1247 for (Index za_inc_index = 0; za_inc_index < nza; ++za_inc_index) {
1248 Complex s11;
1249 Complex s12;
1250 Complex s21;
1251 Complex s22;
1252 VectorView K =
1253 ssd.ext_mat_data(f_index, T_index, za_inc_index, 0, joker);
1254
1255 const Numeric beta = 0.;
1256 const Numeric alpha = 0.;
1257 ampl_(nmax,
1258 lam_f,
1259 ssd.za_grid[za_inc_index],
1260 ssd.za_grid[za_inc_index],
1261 0.,
1262 0.,
1263 alpha,
1264 beta,
1265 s11,
1266 s12,
1267 s21,
1268 s22);
1269
1270 K[0] = (Complex(0., -1.) * (s11 + s22)).real();
1271 K[1] = (Complex(0., 1.) * (s22 - s11)).real();
1272 K[2] = (s22 - s11).real();
1273
1274 K *= lam_f;
1275 }
1276 }
1277 }
1278
1279 csca_data *= 2. * PI * PI / 32400.;
1280 ssd.abs_vec_data =
1281 ssd.ext_mat_data(joker, joker, joker, joker, Range(0, 2));
1282 ssd.abs_vec_data -= csca_data;
1283
1284 break;
1285 }
1286 default: {
1287 std::ostringstream os;
1288 os << "Only particle types totally_random and azimuthally_random\n"
1289 << "are currently supported: " << ssd.ptype;
1290 throw std::runtime_error(os.str());
1291 break;
1292 }
1293 }
1294}
1295
1296// Documentation in header file.
1297void tmatrix_ampld_test(const Verbosity& verbosity) {
1299
1300 out0 << "======================================================\n";
1301 out0 << "Test for nonspherical particles in a fixed orientation\n";
1302 out0 << "Output should match 3rdparty/tmatrix/tmatrix_ampld.ref\n";
1303 out0 << "======================================================\n";
1304
1305 // Same inputs as in example included in original ampld.lp.f
1306 Numeric rat = 1.;
1307 Numeric axi = 10.; //[um]
1308 Index np = -1;
1309 Numeric lam = acos(-1.) * 2.; //[um]
1310 Numeric eps = 0.5;
1311 Numeric mrr = 1.5;
1312 Numeric mri = 0.02;
1313 Numeric ddelt = 0.001;
1314
1315 Index quiet = 1;
1316
1317 // Output variables
1318 Index nmax;
1319 Numeric csca;
1320 Numeric cext;
1321 char errmsg[1024] = "";
1322
1323 tmatrix_(
1324 rat, axi, np, lam, eps, mrr, mri, ddelt, quiet, nmax, csca, cext, errmsg);
1325
1326 out0 << "nmax: " << nmax << "\n";
1327 out0 << "csca: " << csca << " um2\n";
1328 out0 << "cext: " << cext << " um2\n";
1329
1330 out0 << "Error message: " << (strlen(errmsg) ? errmsg : "None") << "\n";
1331
1332 // Same inputs as in example included in original ampld.lp.f
1333 Numeric alpha = 145.;
1334 Numeric beta = 52.;
1335 Numeric thet0 = 56.;
1336 Numeric thet = 65.;
1337 Numeric phi0 = 114.;
1338 Numeric phi = 128.;
1339
1340 // Output variables
1341 Complex s11;
1342 Complex s12;
1343 Complex s21;
1344 Complex s22;
1345 ampl_(nmax, lam, thet0, thet, phi0, phi, alpha, beta, s11, s12, s21, s22);
1346
1347 out0 << "AMPLITUDE MATRIX (all in [um]): \n";
1348 out0 << "s11: " << s11 << "\n";
1349 out0 << "s12: " << s12 << "\n";
1350 out0 << "s21: " << s21 << "\n";
1351 out0 << "s22: " << s22 << "\n";
1352
1353 Matrix z;
1354 ampmat_to_phamat(z, s11, s12, s21, s22);
1355 //z *= 1e12; // meter^2 to micron^2 for comparison with original results
1356
1357 out0 << "PHASE MATRIX (all un [um2]): \n" << z << "\n";
1358}
1359
1360// Documentation in header file.
1361void tmatrix_tmd_test(const Verbosity& verbosity) {
1363
1364 out0 << "======================================================\n";
1365 out0 << "Test for randomly oriented nonspherical particles\n";
1366 out0 << "Output should match 3rdparty/tmatrix/tmatrix_tmd.ref\n";
1367 out0 << "======================================================\n";
1368
1369 // Same inputs as in example included in original tmd.lp.f
1370 Numeric rat = 0.5;
1371 Index ndistr = 3;
1372 Numeric axmax = 1.; //[um]
1373 Index npnax = 2;
1374 Numeric b = 0.1;
1375 Numeric gam = 0.5;
1376 Index nkmax = 5;
1377 Numeric eps = 2;
1378 Index np = -1;
1379 Numeric lam = 0.5; //[um]
1380 Numeric mrr = 1.53;
1381 Numeric mri = 0.008;
1382 Numeric ddelt = 0.001;
1383 Index npna = 19;
1384 Index ndgs = 2;
1385 Numeric r1rat = 0.89031;
1386 Numeric r2rat = 1.56538;
1387
1388 Index quiet = 1;
1389
1390 // Output variables
1391 Numeric reff;
1392 Numeric veff;
1393 Numeric cext;
1394 Numeric csca;
1395 Numeric walb;
1396 Numeric asymm;
1397 Vector f11(npna, 0.);
1398 Vector f22(npna, 0.);
1399 Vector f33(npna, 0.);
1400 Vector f44(npna, 0.);
1401 Vector f12(npna, 0.);
1402 Vector f34(npna, 0.);
1403 char errmsg[1024] = "";
1404
1405 tmd_(rat,
1406 ndistr,
1407 axmax,
1408 npnax,
1409 b,
1410 gam,
1411 nkmax,
1412 eps,
1413 np,
1414 lam,
1415 mrr,
1416 mri,
1417 ddelt,
1418 npna,
1419 ndgs,
1420 r1rat,
1421 r2rat,
1422 quiet,
1423 reff,
1424 veff,
1425 cext,
1426 csca,
1427 walb,
1428 asymm,
1429 f11.get_c_array(),
1430 f22.get_c_array(),
1431 f33.get_c_array(),
1432 f44.get_c_array(),
1433 f12.get_c_array(),
1434 f34.get_c_array(),
1435 errmsg);
1436
1437 out0 << "reff: " << reff << " um\n";
1438 out0 << "veff: " << veff << "\n";
1439 out0 << "cext: " << cext << " um2\n";
1440 out0 << "csca: " << csca << " um2\n";
1441 out0 << "walb: " << walb << "\n";
1442 out0 << "asymm: " << asymm << "\n";
1443 out0 << "f11: " << f11 << "\n";
1444 out0 << "f22: " << f22 << "\n";
1445 out0 << "f33: " << f33 << "\n";
1446 out0 << "f44: " << f44 << "\n";
1447 out0 << "f12: " << f12 << "\n";
1448 out0 << "f34: " << f34 << "\n";
1449
1450 out0 << "Error message: " << (strlen(errmsg) ? errmsg : "None") << "\n";
1451}
1452
1453// Documentation in header file.
1454void calc_ssp_random_test(const Verbosity& verbosity) {
1456 out0 << "======================================================\n";
1457 out0 << "Test calculation of single scattering data\n";
1458 out0 << "for randomly oriented, oblate particles\n";
1459 out0 << "======================================================\n";
1460
1462
1463 ssd.ptype = PTYPE_TOTAL_RND;
1464 ssd.f_grid = {230e9, 240e9};
1465 ssd.T_grid = {220, 250};
1466 nlinspace(ssd.za_grid, 0, 180, 19);
1467 nlinspace(ssd.aa_grid, 0, 180, 19);
1468
1469 // Refractive index real and imagenary parts
1470 // Dimensions: [nf, nT];
1471 Matrix mrr(ssd.f_grid.nelem(), ssd.T_grid.nelem(), 1.78031135);
1472 Matrix mri(ssd.f_grid.nelem(), ssd.T_grid.nelem(), 0.00278706);
1473
1474 mrr(0, 0) = 1.78031135;
1475 mrr(0, 1) = 1.78150475;
1476 mrr(1, 0) = 1.78037238;
1477 mrr(1, 1) = 1.78147686;
1478
1479 mri(0, 0) = 0.00278706;
1480 mri(0, 1) = 0.00507565;
1481 mri(1, 0) = 0.00287245;
1482 mri(1, 1) = 0.00523012;
1483
1484 calcSingleScatteringDataProperties(ssd, mrr, mri, 200.e-6, -1, 1.5);
1485
1486 out0 << "ssd.pha_mat_data(0, 0, joker, 0, 0, joker, joker):\n"
1487 << ssd.pha_mat_data(0, 0, joker, 0, 0, joker, joker) << "\n\n";
1488
1489 out0 << "ssd.ext_mat_data:\n" << ssd.ext_mat_data << "\n\n";
1490 out0 << "ssd.abs_vec_data:\n" << ssd.abs_vec_data << "\n\n";
1491
1492 out0 << "======================================================\n";
1493 out0 << "Test calculation of single scattering data\n";
1494 out0 << "for randomly oriented, prolate particles\n";
1495 out0 << "======================================================\n";
1496
1497 calcSingleScatteringDataProperties(ssd, mrr, mri, 200.e-6, -1, 0.7);
1498
1499 out0 << "ssd.pha_mat_data(0, 0, joker, 0, 0, joker, joker):\n"
1500 << ssd.pha_mat_data(0, 0, joker, 0, 0, joker, joker) << "\n\n";
1501
1502 out0 << "ssd.ext_mat_data:\n" << ssd.ext_mat_data << "\n\n";
1503 out0 << "ssd.abs_vec_data:\n" << ssd.abs_vec_data << "\n\n";
1504}
1505
1506// Documentation in header file.
1507void calc_ssp_fixed_test(const Verbosity& verbosity) {
1509 out0 << "======================================================\n";
1510 out0 << "Test calculation of single scattering data\n";
1511 out0 << "for oblate particles with fixed orientation\n";
1512 out0 << "======================================================\n";
1513
1515
1517 ssd.f_grid = {230e9, 240e9};
1518 ssd.T_grid = {220, 250};
1519 nlinspace(ssd.za_grid, 0, 180, 19);
1520 nlinspace(ssd.aa_grid, 0, 180, 19);
1521
1522 // Refractive index real and imagenary parts
1523 // Dimensions: [nf, nT];
1524 Matrix mrr(ssd.f_grid.nelem(), ssd.T_grid.nelem(), 1.78031135);
1525 Matrix mri(ssd.f_grid.nelem(), ssd.T_grid.nelem(), 0.00278706);
1526
1527 mrr(0, 0) = 1.78031135;
1528 mrr(0, 1) = 1.78150475;
1529 mrr(1, 0) = 1.78037238;
1530 mrr(1, 1) = 1.78147686;
1531
1532 mri(0, 0) = 0.00278706;
1533 mri(0, 1) = 0.00507565;
1534 mri(1, 0) = 0.00287245;
1535 mri(1, 1) = 0.00523012;
1536
1537 calcSingleScatteringDataProperties(ssd, mrr, mri, 200.e-6, -1, 1.5);
1538
1539 out0 << "ssd.pha_mat_data(0, 0, joker, 0, 0, joker, joker):\n"
1540 << ssd.pha_mat_data(0, 0, joker, 0, 0, joker, joker) << "\n\n";
1541
1542 out0 << "ssd.ext_mat_data(0, 0, joker, joker, joker):\n"
1543 << ssd.ext_mat_data(0, 0, joker, joker, joker) << "\n\n";
1544
1545 out0 << "abs_vec_data:\n" << ssd.abs_vec_data << "\n\n";
1546
1547 out0 << "======================================================\n";
1548 out0 << "Test calculation of single scattering data\n";
1549 out0 << "for prolate particles with fixed orientation\n";
1550 out0 << "======================================================\n";
1551 calcSingleScatteringDataProperties(ssd, mrr, mri, 200.e-6, -1, 0.7);
1552
1553 out0 << "ssd.pha_mat_data(0, 0, joker, 0, 0, joker, joker):\n"
1554 << ssd.pha_mat_data(0, 0, joker, 0, 0, joker, joker) << "\n\n";
1555
1556 out0 << "ssd.ext_mat_data(0, 0, joker, joker, joker):\n"
1557 << ssd.ext_mat_data(0, 0, joker, joker, joker) << "\n\n";
1558
1559 out0 << "abs_vec_data:\n" << ssd.abs_vec_data << "\n\n";
1560}
const Numeric SPEED_OF_LIGHT
A constant view of a Matrix.
Definition: matpackI.h:1014
Index nrows() const ARTS_NOEXCEPT
Returns the number of rows.
Definition: matpackI.cc:433
Index ncols() const ARTS_NOEXCEPT
Returns the number of columns.
Definition: matpackI.cc:436
Index nelem() const ARTS_NOEXCEPT
Returns the number of elements.
Definition: matpackI.cc:48
The Matrix class.
Definition: matpackI.h:1225
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
The range class.
Definition: matpackI.h:165
The Tensor5 class.
Definition: matpackV.h:506
void resize(Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackV.cc:1741
void resize(Index l, Index v, Index s, Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackVII.cc:5482
The VectorView class.
Definition: matpackI.h:626
const Numeric * get_c_array() const ARTS_NOEXCEPT
Conversion to plain C-array, const-version.
Definition: matpackI.cc:271
The Vector class.
Definition: matpackI.h:876
void resize(Index n)
Resize function.
Definition: matpackI.cc:408
const Numeric PI
Global constant, pi.
#define beta
#define precision
Definition: logic.cc:43
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
Definition: math_funcs.cc:225
Implementation of Matrix, Vector, and such stuff.
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
const Joker joker
constexpr Complex conj(Complex c) noexcept
the conjugate of c
std::complex< Numeric > Complex
constexpr Numeric real(Complex c) noexcept
real
Declarations having to do with the four output streams.
#define CREATE_OUT0
Definition: messages.h:204
Scattering database structure and functions.
@ PTYPE_AZIMUTH_RND
Definition: optproperties.h:54
@ PTYPE_TOTAL_RND
Definition: optproperties.h:55
#define c
#define b
void calc_ssp_fixed_test(const Verbosity &verbosity)
Single scattering properties calculation for particles with fixed orientation.
Definition: tmatrix.cc:1507
void calcSingleScatteringDataProperties(SingleScatteringData &ssd, ConstMatrixView ref_index_real, ConstMatrixView ref_index_imag, const Numeric equiv_radius, const Index np, const Numeric aspect_ratio, const Numeric precision, const Index ndgs, const Index robust, const Index quiet)
Calculate SingleScatteringData properties.
Definition: tmatrix.cc:979
void integrate_phamat_alpha6(Matrix &phamat, const Index &nmax, const Numeric &lam, const Numeric &thet0, const Numeric &thet, const Numeric &phi0, const Numeric &phi, const Numeric &beta, const Numeric &alpha_1, const Numeric &alpha_2)
Integrate phase matrix over particle orientation angle.
Definition: tmatrix.cc:560
void tmatrix_fixed_orientation(Numeric &cext, Numeric &csca, Index &nmax, const Numeric equiv_radius, const Numeric aspect_ratio, const Index np, const Numeric lam, const Numeric ref_index_real, const Numeric ref_index_imag, const Numeric precision, const Index quiet=1)
Calculate properties for particles in a fixed orientation.
Definition: tmatrix.cc:941
void calc_phamat(Matrix &z, const Index &nmax, const Numeric &lam, const Numeric &thet0, const Numeric &thet, const Numeric &phi0, const Numeric &phi, const Numeric &beta, const Numeric &alpha)
Definition: tmatrix.cc:419
void tmatrix_random_orientation(Numeric &cext, Numeric &csca, Vector &f11, Vector &f22, Vector &f33, Vector &f44, Vector &f12, Vector &f34, const Numeric equiv_radius, const Numeric aspect_ratio, const Index np, const Numeric lam, const Numeric ref_index_real, const Numeric ref_index_imag, const Numeric precision, const Index nza, const Index ndgs, const Index quiet=1)
Calculate properties for randomly oriented particles.
Definition: tmatrix.cc:843
void tmatrix_tmd_test(const Verbosity &verbosity)
T-Matrix validation test.
Definition: tmatrix.cc:1361
void tmatrix_(const Numeric &rat, const Numeric &axi, const Index &np, const Numeric &lam, const Numeric &eps, const Numeric &mrr, const Numeric &mri, const Numeric &ddelt, const Index &quiet, Index &nmax, Numeric &csca, Numeric &cext, char *errmsg)
T-matrix code for nonspherical particles in a fixed orientation.
Definition: tmatrix.cc:378
void calc_ssp_random_test(const Verbosity &verbosity)
Single scattering properties calculation for randomly oriented particles.
Definition: tmatrix.cc:1454
void avgtmatrix_(const Index &nmax)
Perform orientation averaging.
Definition: tmatrix.cc:412
void tmatrix_ampld_test(const Verbosity &verbosity)
T-Matrix validation test.
Definition: tmatrix.cc:1297
void ampl_(const Index &nmax, const Numeric &lam, const Numeric &thet0, const Numeric &thet, const Numeric &phi0, const Numeric &phi, const Numeric &alpha, const Numeric &beta, Complex &s11, Complex &s12, Complex &s21, Complex &s22)
T-matrix code for nonspherical particles in a fixed orientation.
Definition: tmatrix.cc:396
void tmd_(const Numeric &rat, const Index &ndistr, const Numeric &axmax, const Index &npnax, const Numeric &b, const Numeric &gam, const Index &nkmax, const Numeric &eps, const Index &np, const Numeric &lam, const Numeric &mrr, const Numeric &mri, const Numeric &ddelt, const Index &npna, const Index &ndgs, const Numeric &r1rat, const Numeric &r2rat, const Index &quiet, Numeric &reff, Numeric &veff, Numeric &cext, Numeric &csca, Numeric &walb, Numeric &asymm, Numeric *f11, Numeric *f22, Numeric *f33, Numeric *f44, Numeric *f12, Numeric *f34, char *errmsg)
T-matrix code for randomly oriented nonspherical particles.
Definition: tmatrix.cc:342
void integrate_phamat_alpha10(Matrix &phamat, const Index &nmax, const Numeric &lam, const Numeric &thet0, const Numeric &thet, const Numeric &phi0, const Numeric &phi, const Numeric &beta, const Numeric &alpha_1, const Numeric &alpha_2)
Integrate phase matrix over particle orientation angle.
Definition: tmatrix.cc:510
void ampmat_to_phamat(Matrix &z, const Complex &s11, const Complex &s12, const Complex &s21, const Complex &s22)
Calculate phase matrix from amplitude matrix.
Definition: tmatrix.cc:450
void integrate_phamat_theta0_phi10(Matrix &phamat, const Index &nmax, const Numeric &lam, const Numeric &thet0_1, const Numeric &thet0_2, const Numeric &thet, const Numeric &phi0, const Numeric &phi_1, const Numeric &phi_2, const Numeric &beta, const Numeric &alpha)
Integrate phase matrix over angles thet0 and phi.
Definition: tmatrix.cc:612
void integrate_phamat_theta0_phi_alpha6(Matrix &phamat, const Index &nmax, const Numeric &lam, const Numeric &thet0_1, const Numeric &thet0_2, const Numeric &thet, const Numeric &phi0, const Numeric &phi_1, const Numeric &phi_2, const Numeric &beta, const Numeric &alpha_1, const Numeric &alpha_2)
Integrate phase matrix over angles thet0, phi and alpha.
Definition: tmatrix.cc:723
Declarations for the T-Matrix interface.