ARTS 2.5.11 (git: 6827797f)
transmissionmatrix.h
Go to the documentation of this file.
1
11#ifndef transmissionmatrix_h
12#define transmissionmatrix_h
13
14#pragma GCC diagnostic push
15#pragma GCC diagnostic ignored "-Wshadow"
16#pragma GCC diagnostic ignored "-Wconversion"
17#pragma GCC diagnostic ignored "-Wfloat-conversion"
18#if defined(__clang__)
19#pragma GCC diagnostic ignored "-Wdeprecated-copy-with-dtor"
20#pragma GCC diagnostic ignored "-Wdeprecated-copy-with-user-provided-dtor"
21#else
22#pragma GCC diagnostic ignored "-Wclass-memaccess"
23#endif
24
25#include <Eigen/Dense>
26
27#pragma GCC diagnostic pop
28
29#include "jacobian.h"
30#include "propagationmatrix.h"
31
32template <int N>
33Eigen::Matrix<Numeric, N, N> prop_matrix(const ConstVectorView& vec) {
34#define a vec[0]
35#define b vec[1]
36#define c vec[2]
37#define d vec[3]
38#define u vec[N]
39#define v vec[5]
40#define w vec[6]
41 static_assert (N>0 and N<5, "Bad size N");
42 if constexpr (N == 1) {
43 return Eigen::Matrix<Numeric, 1, 1>(a);
44 } else if constexpr (N==2) {
45 return (Eigen::Matrix2d() << a, b, b, a).finished();
46 } else if constexpr (N==3) {
47 return (Eigen::Matrix3d() << a, b, c, b, a, u, c, -u, a).finished();
48 } else if constexpr (N==4) {
49 return (Eigen::Matrix4d() << a, b, c, d, b, a, u, v, c, -u, a, w, d, -v, -w, a).finished();
50 }
51#undef a
52#undef b
53#undef c
54#undef d
55#undef u
56#undef v
57#undef w
58}
59
60template <int N>
61Eigen::Matrix<Numeric, N, N> prop_matrix(const ConstMatrixView& m) {
62 static_assert (N>0 and N<5, "Bad size N");
63 if constexpr (N == 1) {
64 return Eigen::Matrix<Numeric, 1, 1>(m(0, 0));
65 } else if constexpr (N == 2) {
66 return (Eigen::Matrix2d() << m(0, 0), m(0, 1),
67 m(1, 0), m(1, 1)).finished();
68 } else if constexpr (N == 3) {
69 return (Eigen::Matrix3d() << m(0, 0), m(0, 1), m(0, 2),
70 m(1, 0), m(1, 1), m(1, 2),
71 m(2, 0), m(2, 1), m(2, 2)) .finished();
72 } else if constexpr (N == 4) {
73 return (Eigen::Matrix4d() << m(0, 0), m(0, 1), m(0, 2), m(0, 3),
74 m(1, 0), m(1, 1), m(1, 2), m(1, 3),
75 m(2, 0), m(2, 1), m(2, 2), m(2, 3),
76 m(3, 0), m(3, 1), m(3, 2), m(3, 3)) .finished();
77 }
78}
79
80template <int N>
81Eigen::Matrix<Numeric, N, N> inv_prop_matrix(const ConstVectorView& vec) {
82#define a vec[0]
83#define b vec[1]
84#define c vec[2]
85#define d vec[3]
86#define u vec[N]
87#define v vec[5]
88#define w vec[6]
89 static_assert (N>0 and N<5, "Bad size N");
90 if constexpr (N == 1) {
91 return Eigen::Matrix<double, 1, 1>(1 / a);
92 } else if constexpr (N==2) {
93 return (Eigen::Matrix2d() << a, -b, -b, a).finished() / (a * a - b * b);
94 } else if constexpr (N==3) {
95 return (Eigen::Matrix3d() << a * a + u * u,
96 -a * b - c * u,
97 -a * c + b * u,
98 -a * b + c * u,
99 a * a - c * c,
100 -a * u + b * c,
101 -a * c - b * u,
102 a * u + b * c,
103 a * a - b * b)
104 .finished() /
105 (a * a * a - a * b * b - a * c * c + a * u * u);
106 } else if constexpr (N==4) {
107 return (Eigen::Matrix4d() << a * a * a + a * u * u + a * v * v + a * w * w,
108 -a * a * b - a * c * u - a * d * v - b * w * w + c * v * w -
109 d * u * w,
110 -a * a * c + a * b * u - a * d * w + b * v * w - c * v * v +
111 d * u * v,
112 -a * a * d + a * b * v + a * c * w - b * u * w + c * u * v -
113 d * u * u,
114 -a * a * b + a * c * u + a * d * v - b * w * w + c * v * w -
115 d * u * w,
116 a * a * a - a * c * c - a * d * d + a * w * w,
117 -a * a * u + a * b * c - a * v * w + b * d * w - c * d * v +
118 d * d * u,
119 -a * a * v + a * b * d + a * u * w - b * c * w + c * c * v -
120 c * d * u,
121 -a * a * c - a * b * u + a * d * w + b * v * w - c * v * v +
122 d * u * v,
123 a * a * u + a * b * c - a * v * w - b * d * w + c * d * v - d * d * u,
124 a * a * a - a * b * b - a * d * d + a * v * v,
125 -a * a * w + a * c * d - a * u * v + b * b * w - b * c * v +
126 b * d * u,
127 -a * a * d - a * b * v - a * c * w - b * u * w + c * u * v -
128 d * u * u,
129 a * a * v + a * b * d + a * u * w + b * c * w - c * c * v + c * d * u,
130 a * a * w + a * c * d - a * u * v - b * b * w + b * c * v - b * d * u,
131 a * a * a - a * b * b - a * c * c + a * u * u)
132 .finished() /
133 (a * a * a * a - a * a * b * b - a * a * c * c - a * a * d * d +
134 a * a * u * u + a * a * v * v + a * a * w * w - b * b * w * w +
135 2 * b * c * v * w - 2 * b * d * u * w - c * c * v * v +
136 2 * c * d * u * v - d * d * u * u);
137 }
138#undef a
139#undef b
140#undef c
141#undef d
142#undef u
143#undef v
144#undef w
145}
146
150 std::vector<Eigen::Matrix4d> T4;
151 std::vector<Eigen::Matrix3d> T3;
152 std::vector<Eigen::Matrix2d> T2;
153 std::vector<Eigen::Matrix<double, 1, 1>> T1;
154
160 TransmissionMatrix(Index nf = 0, Index stokes = 1);
161
166 TransmissionMatrix(TransmissionMatrix&& tm) noexcept = default;
167
173
180
187
193 explicit TransmissionMatrix(const PropagationMatrix& pm,
194 const Numeric& r = 1.0);
195
196 operator Tensor3() const;
197
198 explicit TransmissionMatrix(const ConstMatrixView& mat);
199
205 [[nodiscard]] const Eigen::Matrix4d& Mat4(size_t i) const;
206
212 [[nodiscard]] const Eigen::Matrix3d& Mat3(size_t i) const;
213
219 [[nodiscard]] const Eigen::Matrix2d& Mat2(size_t i) const;
220
226 [[nodiscard]] const Eigen::Matrix<double, 1, 1>& Mat1(size_t i) const;
227
233 [[nodiscard]] Eigen::MatrixXd Mat(size_t i) const;
234
240 Eigen::Matrix4d& Mat4(size_t i);
241
247 Eigen::Matrix3d& Mat3(size_t i);
248
254 Eigen::Matrix2d& Mat2(size_t i);
255
261 Eigen::Matrix<double, 1, 1>& Mat1(size_t i);
262
264 void setIdentity();
265
267 void setZero();
268
276 void mul(const TransmissionMatrix& A, const TransmissionMatrix& B);
277
285 void mul_aliased(const TransmissionMatrix& A, const TransmissionMatrix& B);
286
294 [[nodiscard]] Numeric operator()(const Index i,
295 const Index j,
296 const Index k) const;
297
305 Numeric& operator()(const Index i, const Index j, const Index k);
306
308 [[nodiscard]] Index Frequencies() const;
309
315 TransmissionMatrix& operator+=(const LazyScale<TransmissionMatrix>& lstm);
316
322 TransmissionMatrix& operator*=(const Numeric& scale);
323
329 TransmissionMatrix& operator=(const LazyScale<TransmissionMatrix>& lstm);
330
332 friend std::ostream& operator<<(std::ostream& os,
333 const TransmissionMatrix& tm);
334
336 friend std::istream& operator>>(std::istream& data, TransmissionMatrix& tm);
337
339 template <int N>
340 auto& TraMat(size_t i) noexcept {
341 static_assert(N > 0 and N < 5, "Bad size N");
342 if constexpr (N == 1)
343 return T1[i];
344 else if constexpr (N == 2)
345 return T2[i];
346 else if constexpr (N == 3)
347 return T3[i];
348 else if constexpr (N == 4)
349 return T4[i];
350 }
351
353 template <int N>
354 [[nodiscard]] auto& TraMat(size_t i) const noexcept {
355 static_assert(N > 0 and N < 5, "Bad size N");
356 if constexpr (N == 1)
357 return T1[i];
358 else if constexpr (N == 2)
359 return T2[i];
360 else if constexpr (N == 3)
361 return T3[i];
362 else if constexpr (N == 4)
363 return T4[i];
364 }
365
367 template <int N>
368 [[nodiscard]] Eigen::Matrix<Numeric, N, N> OptDepth(size_t i) const noexcept {
369 return (-TraMat<N>(i).diagonal().array().log().matrix()).asDiagonal();
370 }
371
391 template <int N>
392 [[nodiscard]] Eigen::Matrix<Numeric, N, 1> second_order_integration_source(
393 const Eigen::Matrix<Numeric, N, N> T,
394 const Eigen::Matrix<Numeric, N, 1> far,
395 const Eigen::Matrix<Numeric, N, 1> close,
396 const Eigen::Matrix<Numeric, N, N> Kfar,
397 const Eigen::Matrix<Numeric, N, N> Kclose,
398 const Numeric r) const noexcept {
399 static_assert(N > 0 and N < 5, "Bad size N");
400
401 const auto I = Eigen::Matrix<Numeric, N, N>::Identity();
402 if (T(0, 0) < 0.99) {
403 const auto od = 0.5 * r * (Kfar + Kclose);
404 Eigen::Matrix<Numeric, N, 1> second =
405 od.inverse() * ((I - (I + od) * T) * far + (od - I + T) * close);
406 if ((far[0] < close[0] and Kfar(0, 0) > Kclose(0, 0)) or
407 (far[0] > close[0] and Kfar(0, 0) < Kclose(0, 0))) {
408 Eigen::Matrix<Numeric, N, 1> second_limit =
409 0.5 * (I - T) * (far + close);
410
411 // FIXME: This is the equation given in the source material...
412 // const Eigen::Matrix<Numeric, N, 1> second_limit = 0.25 * r * (Kfar * far + Kclose * close);
413
414 if (second_limit[0] > second[0]) return second;
415 return second_limit;
416 }
417
418 return second;
419 }
420
421 return 0.5 * (I - T) * (far + close);
422 }
423
424 template <int N>
425 [[nodiscard]] Eigen::Matrix<Numeric, N, 1> second_order_integration_dsource(
426 size_t i,
427 const TransmissionMatrix& dx,
428 const Eigen::Matrix<Numeric, N, 1> far,
429 const Eigen::Matrix<Numeric, N, 1> close,
430 const Eigen::Matrix<Numeric, N, 1> d,
431 bool isfar) const noexcept {
432 static_assert(N > 0 and N < 5, "Bad size N");
433
434 const auto I = Eigen::Matrix<Numeric, N, N>::Identity();
435 const auto T = TraMat<N>(i);
436 const auto& dTdx = dx.TraMat<N>(i);
437 const Eigen::Matrix<Numeric, N, 1> first = 0.5 * (I - T) * (far + close);
438 if (T(0, 0) < 0.99) {
439 const auto od = OptDepth<N>(i);
440 const auto doddx = dx.OptDepth<N>(i);
441 const Eigen::Matrix<Numeric, N, 1> second =
442 od.inverse() * ((I - (I + od) * T) * far + (od - I + T) * close);
443 if (first[0] > second[0]) {
444 if (isfar) {
445 return od.inverse() *
446 ((I - (I + od) * T) * d - (I + od) * dTdx * far +
447 (doddx + T) * close - doddx * second);
448 }
449
450 return od.inverse() * (-(I + od) * dTdx * far + (od - I + T) * d +
451 (doddx + T) * close - doddx * second);
452 }
453
454 return 0.5 * ((I - T) * d - dTdx * (far + close));
455 }
456
457 return 0.5 * ((I - T) * d - dTdx * (far + close));
458 }
459};
460
467[[nodiscard]] LazyScale<TransmissionMatrix> operator*(
468 const TransmissionMatrix& tm, const Numeric& x);
469
476[[nodiscard]] LazyScale<TransmissionMatrix> operator*(
477 const Numeric& x, const TransmissionMatrix& tm);
478
482 std::vector<Eigen::Vector4d> R4;
483 std::vector<Eigen::Vector3d> R3;
484 std::vector<Eigen::Vector2d> R2;
485 std::vector<Eigen::Matrix<double, 1, 1>> R1;
486
492 RadiationVector(Index nf = 0, Index stokes = 1);
493
498 RadiationVector(RadiationVector&& rv) noexcept = default;
499
504 RadiationVector(const RadiationVector& rv) = default;
505
512
518 RadiationVector& operator=(RadiationVector&& rv) noexcept = default;
519
526
531 void leftMul(const TransmissionMatrix& T);
532
537 void SetZero(size_t i);
538
539 void SetZero();
540
546 [[nodiscard]] const Eigen::Vector4d& Vec4(size_t i) const;
547
553 [[nodiscard]] const Eigen::Vector3d& Vec3(size_t i) const;
554
560 [[nodiscard]] const Eigen::Vector2d& Vec2(size_t i) const;
561
567 [[nodiscard]] const Eigen::Matrix<double, 1, 1>& Vec1(size_t i) const;
568
574 [[nodiscard]] Eigen::VectorXd Vec(size_t i) const;
575
581 Eigen::Vector4d& Vec4(size_t i);
582
588 Eigen::Vector3d& Vec3(size_t i);
589
595 Eigen::Vector2d& Vec2(size_t i);
596
602 Eigen::Matrix<double, 1, 1>& Vec1(size_t i);
603
609 void rem_avg(const RadiationVector& O1, const RadiationVector& O2);
610
616 void add_avg(const RadiationVector& O1, const RadiationVector& O2);
617
624 void add_weighted(const TransmissionMatrix& T,
625 const RadiationVector& far,
626 const RadiationVector& close,
627 const ConstMatrixView& Kfar,
628 const ConstMatrixView& Kclose,
629 const Numeric r);
630
639 void addDerivEmission(const TransmissionMatrix& PiT,
640 const TransmissionMatrix& dT,
641 const TransmissionMatrix& T,
642 const RadiationVector& ImJ,
643 const RadiationVector& dJ);
644
654 const TransmissionMatrix& dT,
655 const TransmissionMatrix& T,
656 const RadiationVector& I,
657 const RadiationVector& far,
658 const RadiationVector& close,
659 const RadiationVector& d,
660 bool isfar);
661
669 const TransmissionMatrix& dT,
670 const RadiationVector& I);
671
681 void addMultiplied(const TransmissionMatrix& A, const RadiationVector& x);
682
691 const TransmissionMatrix& PiT,
692 const TransmissionMatrix& Z,
693 const TransmissionMatrix& dZ);
694
703 const TransmissionMatrix& Tr,
704 const TransmissionMatrix& Tf,
705 const TransmissionMatrix& Z);
706
715 const TransmissionMatrix& Tr,
716 const TransmissionMatrix& Tf,
717 const TransmissionMatrix& dZ);
718
724 RadiationVector& operator=(const ConstMatrixView& M);
725
732 const Numeric& operator()(const Index i, const Index j) const;
733
740 [[nodiscard]] Numeric& operator()(const Index i, const Index j);
745 operator Matrix() const;
746
754 void setSource(const StokesVector& a,
755 const ConstVectorView& B,
756 const StokesVector& S,
757 Index i);
758
760 [[nodiscard]] Index Frequencies() const;
761
763 friend std::ostream& operator<<(std::ostream& os, const RadiationVector& rv);
764
766 friend std::istream& operator>>(std::istream& data, RadiationVector& rv);
767};
768
777
782};
783
786 Forward,
787 Reverse,
788};
789
792 Emission,
795};
796
815 const RadiationVector& J1,
816 const RadiationVector& J2,
817 const ArrayOfRadiationVector& dJ1,
818 const ArrayOfRadiationVector& dJ2,
819 const TransmissionMatrix& T,
820 const TransmissionMatrix& PiT,
821 const ArrayOfTransmissionMatrix& dT1,
822 const ArrayOfTransmissionMatrix& dT2,
823 const PropagationMatrix& K1,
824 const PropagationMatrix& K2,
825 const ArrayOfPropagationMatrix& dK1,
826 const ArrayOfPropagationMatrix& dK2,
827 const Numeric r,
828 const Vector& dr1,
829 const Vector& dr2,
830 const Index ia,
831 const Index iz,
832 const RadiativeTransferSolver solver);
833
852 RadiationVector& J_add,
853 const PropagationMatrix& K,
854 const StokesVector& a,
855 const StokesVector& S,
856 const ArrayOfPropagationMatrix& dK,
857 const ArrayOfStokesVector& da,
858 const ArrayOfStokesVector& dS,
859 const ConstVectorView& B,
860 const ConstVectorView& dB_dT,
861 const ArrayOfRetrievalQuantity& jacobian_quantities,
862 const bool& jacobian_do);
863
881 const PropagationMatrix& K1,
882 const PropagationMatrix& K2,
883 const ArrayOfPropagationMatrix& dK1,
884 const ArrayOfPropagationMatrix& dK2,
885 const Numeric& r,
886 const Numeric& dr_dtemp1,
887 const Numeric& dr_dtemp2,
888 const Index temp_deriv_pos);
889
898 const CumulativeTransmission type) /*[[expects: T.nelem()>0]]*/;
899
917 const RadiationVector& I_incoming,
919 const ArrayOfTransmissionMatrix& PiTf,
920 const ArrayOfTransmissionMatrix& PiTr,
925 const BackscatterSolver solver);
926
938ArrayOfTransmissionMatrix bulk_backscatter(const ConstTensor5View& Pe,
939 const ConstMatrixView& pnd);
940
951 const ConstTensor5View& Pe, const ArrayOfMatrix& dpnd_dx);
952
953#endif // transmissionmatrix_h
This can be used to make arrays out of anything.
Definition: array.h:31
Routines for setting up the jacobian.
Radiation Vector for Stokes dimension 1-4.
std::vector< Eigen::Vector2d > R2
Index Frequencies() const
Get frequency count.
RadiationVector(RadiationVector &&rv) noexcept=default
Construct a new Radiation Vector object.
const Eigen::Matrix< double, 1, 1 > & Vec1(size_t i) const
Return Vector at position.
void setSource(const StokesVector &a, const ConstVectorView &B, const StokesVector &S, Index i)
Set this to source vector at position.
void addMultiplied(const TransmissionMatrix &A, const RadiationVector &x)
Add multiply.
const Eigen::Vector4d & Vec4(size_t i) const
Return Vector at position.
void leftMul(const TransmissionMatrix &T)
Multiply radiation vector from the left.
const Numeric & operator()(const Index i, const Index j) const
Access operator.
void rem_avg(const RadiationVector &O1, const RadiationVector &O2)
Remove the average of two other RadiationVector from *this.
std::vector< Eigen::Matrix< double, 1, 1 > > R1
void add_avg(const RadiationVector &O1, const RadiationVector &O2)
Add the average of two other RadiationVector to *this.
std::vector< Eigen::Vector4d > R4
RadiationVector & operator+=(const RadiationVector &rv)
Addition operator.
Eigen::VectorXd Vec(size_t i) const
Return Vector at position by copy.
void addDerivTransmission(const TransmissionMatrix &PiT, const TransmissionMatrix &dT, const RadiationVector &I)
Add the transmission derivative to this.
friend std::ostream & operator<<(std::ostream &os, const RadiationVector &rv)
Output operator.
void setDerivReflection(const RadiationVector &I, const TransmissionMatrix &PiT, const TransmissionMatrix &Z, const TransmissionMatrix &dZ)
Sets *this to the reflection derivative.
const Eigen::Vector3d & Vec3(size_t i) const
Return Vector at position.
const Eigen::Vector2d & Vec2(size_t i) const
Return Vector at position.
std::vector< Eigen::Vector3d > R3
void addWeightedDerivEmission(const TransmissionMatrix &PiT, const TransmissionMatrix &dT, const TransmissionMatrix &T, const RadiationVector &I, const RadiationVector &far, const RadiationVector &close, const RadiationVector &d, bool isfar)
Add the emission derivative to this.
void setBackscatterTransmission(const RadiationVector &I0, const TransmissionMatrix &Tr, const TransmissionMatrix &Tf, const TransmissionMatrix &Z)
Set this to backscatter transmission.
void addDerivEmission(const TransmissionMatrix &PiT, const TransmissionMatrix &dT, const TransmissionMatrix &T, const RadiationVector &ImJ, const RadiationVector &dJ)
Add the emission derivative to this.
RadiationVector & operator=(const RadiationVector &rv)=default
Assign old radiation vector to this.
void add_weighted(const TransmissionMatrix &T, const RadiationVector &far, const RadiationVector &close, const ConstMatrixView &Kfar, const ConstMatrixView &Kclose, const Numeric r)
Add the weighted source of two RadiationVector to *this.
friend std::istream & operator>>(std::istream &data, RadiationVector &rv)
Input operator.
RadiationVector & operator=(RadiationVector &&rv) noexcept=default
Assign old radiation vector to this.
void setBackscatterTransmissionDerivative(const RadiationVector &I0, const TransmissionMatrix &Tr, const TransmissionMatrix &Tf, const TransmissionMatrix &dZ)
Set this to backscatter transmission scatter derivative.
RadiationVector(const RadiationVector &rv)=default
Construct a new Radiation Vector object.
Class to keep track of Transmission Matrices for Stokes Dim 1-4.
TransmissionMatrix & operator=(TransmissionMatrix &&tm) noexcept=default
Move operator.
void mul_aliased(const TransmissionMatrix &A, const TransmissionMatrix &B)
Set this to a multiple of A by B.
std::vector< Eigen::Matrix< double, 1, 1 > > T1
std::vector< Eigen::Matrix4d > T4
Eigen::Matrix< Numeric, N, N > OptDepth(size_t i) const noexcept
Simple template access for the optical depth.
TransmissionMatrix & operator*=(const Numeric &scale)
Scale self.
const Eigen::Matrix2d & Mat2(size_t i) const
Get Matrix at position.
TransmissionMatrix(const TransmissionMatrix &tm)=default
Construct a new Transmission Matrix object.
Eigen::MatrixXd Mat(size_t i) const
Get Matrix at position by copy.
auto & TraMat(size_t i) noexcept
Simple template access for the transmission.
Index Frequencies() const
Number of frequencies.
Eigen::Matrix< Numeric, N, 1 > second_order_integration_dsource(size_t i, const TransmissionMatrix &dx, const Eigen::Matrix< Numeric, N, 1 > far, const Eigen::Matrix< Numeric, N, 1 > close, const Eigen::Matrix< Numeric, N, 1 > d, bool isfar) const noexcept
auto & TraMat(size_t i) const noexcept
Simple template access for the transmission.
void mul(const TransmissionMatrix &A, const TransmissionMatrix &B)
Set this to a multiple of A by B.
std::vector< Eigen::Matrix3d > T3
const Eigen::Matrix4d & Mat4(size_t i) const
Get Matrix at position.
Eigen::Matrix< Numeric, N, 1 > second_order_integration_source(const Eigen::Matrix< Numeric, N, N > T, const Eigen::Matrix< Numeric, N, 1 > far, const Eigen::Matrix< Numeric, N, 1 > close, const Eigen::Matrix< Numeric, N, N > Kfar, const Eigen::Matrix< Numeric, N, N > Kclose, const Numeric r) const noexcept
void setIdentity()
Set to identity matrix.
const Eigen::Matrix< double, 1, 1 > & Mat1(size_t i) const
Get Matrix at position.
std::vector< Eigen::Matrix2d > T2
friend std::istream & operator>>(std::istream &data, TransmissionMatrix &tm)
Input operator.
const Eigen::Matrix3d & Mat3(size_t i) const
Get Matrix at position.
TransmissionMatrix(TransmissionMatrix &&tm) noexcept=default
Construct a new Transmission Matrix object.
friend std::ostream & operator<<(std::ostream &os, const TransmissionMatrix &tm)
Output operator.
Numeric operator()(const Index i, const Index j, const Index k) const
Access value in matrix.
TransmissionMatrix & operator=(const TransmissionMatrix &tm)=default
Assignment operator.
TransmissionMatrix & operator+=(const LazyScale< TransmissionMatrix > &lstm)
Assign to *this lazily.
void setZero()
Set to zero matrix.
RadiativeTransferSolver
Intended to hold various forward solvers.
#define u
ArrayOfTransmissionMatrix bulk_backscatter(const ConstTensor5View &Pe, const ConstMatrixView &pnd)
Bulk back-scattering.
#define d
BackscatterSolver
Intended to hold various backscatter solvers.
CumulativeTransmission
Intended to hold various ways to accumulate the transmission matrix.
#define v
ArrayOfArrayOfTransmissionMatrix bulk_backscatter_derivative(const ConstTensor5View &Pe, const ArrayOfMatrix &dpnd_dx)
Derivatives of bulk back-scattering
LazyScale< TransmissionMatrix > operator*(const TransmissionMatrix &tm, const Numeric &x)
Lazy scale of Transmission Matrix.
ArrayOfTransmissionMatrix cumulative_transmission(const ArrayOfTransmissionMatrix &T, const CumulativeTransmission type)
Accumulate the transmission matrix over all layers.
void set_backscatter_radiation_vector(ArrayOfRadiationVector &I, ArrayOfArrayOfArrayOfRadiationVector &dI, const RadiationVector &I_incoming, const ArrayOfTransmissionMatrix &T, const ArrayOfTransmissionMatrix &PiTf, const ArrayOfTransmissionMatrix &PiTr, const ArrayOfTransmissionMatrix &Z, const ArrayOfArrayOfTransmissionMatrix &dT1, const ArrayOfArrayOfTransmissionMatrix &dT2, const ArrayOfArrayOfTransmissionMatrix &dZ, const BackscatterSolver solver)
Set the backscatter radiation vector.
void stepwise_source(RadiationVector &J, ArrayOfRadiationVector &dJ, RadiationVector &J_add, const PropagationMatrix &K, const StokesVector &a, const StokesVector &S, const ArrayOfPropagationMatrix &dK, const ArrayOfStokesVector &da, const ArrayOfStokesVector &dS, const ConstVectorView &B, const ConstVectorView &dB_dT, const ArrayOfRetrievalQuantity &jacobian_quantities, const bool &jacobian_do)
Set the stepwise source.
#define w
#define a
Eigen::Matrix< Numeric, N, N > inv_prop_matrix(const ConstVectorView &vec)
#define c
#define b
void update_radiation_vector(RadiationVector &I, ArrayOfRadiationVector &dI1, ArrayOfRadiationVector &dI2, const RadiationVector &J1, const RadiationVector &J2, const ArrayOfRadiationVector &dJ1, const ArrayOfRadiationVector &dJ2, const TransmissionMatrix &T, const TransmissionMatrix &PiT, const ArrayOfTransmissionMatrix &dT1, const ArrayOfTransmissionMatrix &dT2, const PropagationMatrix &K1, const PropagationMatrix &K2, const ArrayOfPropagationMatrix &dK1, const ArrayOfPropagationMatrix &dK2, const Numeric r, const Vector &dr1, const Vector &dr2, const Index ia, const Index iz, const RadiativeTransferSolver solver)
Update the Radiation Vector.
void stepwise_transmission(TransmissionMatrix &T, ArrayOfTransmissionMatrix &dT1, ArrayOfTransmissionMatrix &dT2, const PropagationMatrix &K1, const PropagationMatrix &K2, const ArrayOfPropagationMatrix &dK1, const ArrayOfPropagationMatrix &dK2, const Numeric &r, const Numeric &dr_dtemp1, const Numeric &dr_dtemp2, const Index temp_deriv_pos)
Set the stepwise transmission matrix.
Eigen::Matrix< Numeric, N, N > prop_matrix(const ConstVectorView &vec)