29#ifndef transmissionmatrix_h
30#define transmissionmatrix_h
32#pragma GCC diagnostic push
33#pragma GCC diagnostic ignored "-Wshadow"
34#pragma GCC diagnostic ignored "-Wconversion"
35#pragma GCC diagnostic ignored "-Wfloat-conversion"
37#include <Eigen/StdVector>
38#pragma GCC diagnostic pop
46 std::vector<Eigen::Matrix4d, Eigen::aligned_allocator<Eigen::Matrix4d>>
T4;
47 std::vector<Eigen::Matrix3d, Eigen::aligned_allocator<Eigen::Matrix3d>>
T3;
48 std::vector<Eigen::Matrix2d, Eigen::aligned_allocator<Eigen::Matrix2d>>
T2;
49 std::vector<Eigen::Matrix<double, 1, 1>,
50 Eigen::aligned_allocator<Eigen::Matrix<double, 1, 1>>>
76 T1(std::move(tm.T1)) {}
98 T4 = std::move(tm.T4);
99 T3 = std::move(tm.T3);
100 T2 = std::move(tm.T2);
101 T1 = std::move(tm.T1);
114 for (
size_t i = 0; i <
T4.size(); i++)
115 for (
size_t j = 0; j < 4; j++)
116 for (
size_t k = 0; k < 4; k++) T(i, j, k) =
T4[i](j, k);
117 for (
size_t i = 0; i <
T3.size(); i++)
118 for (
size_t j = 0; j < 3; j++)
119 for (
size_t k = 0; k < 3; k++) T(i, j, k) =
T3[i](j, k);
120 for (
size_t i = 0; i <
T2.size(); i++)
121 for (
size_t j = 0; j < 2; j++)
122 for (
size_t k = 0; k < 2; k++) T(i, j, k) =
T2[i](j, k);
123 for (
size_t i = 0; i <
T1.size(); i++) T(i, 0, 0) =
T1[i](0, 0);
132 [[nodiscard]]
const Eigen::Matrix4d&
Mat4(
size_t i)
const {
return T4[i]; }
139 [[nodiscard]]
const Eigen::Matrix3d&
Mat3(
size_t i)
const {
return T3[i]; }
146 [[nodiscard]]
const Eigen::Matrix2d&
Mat2(
size_t i)
const {
return T2[i]; }
153 [[nodiscard]]
const Eigen::Matrix<double, 1, 1>&
Mat1(
size_t i)
const {
return T1[i]; }
160 #pragma GCC diagnostic push
161 #pragma GCC diagnostic ignored "-Wreturn-type"
162 [[nodiscard]] Eigen::MatrixXd
Mat(
size_t i)
const {
174 #pragma GCC diagnostic pop
181 Eigen::Matrix4d&
Mat4(
size_t i) {
return T4[i]; }
188 Eigen::Matrix3d&
Mat3(
size_t i) {
return T3[i]; }
195 Eigen::Matrix2d&
Mat2(
size_t i) {
return T2[i]; }
202 Eigen::Matrix<double, 1, 1>&
Mat1(
size_t i) {
return T1[i]; }
214 std::fill(
T4.begin(),
T4.end(), Eigen::Matrix4d::Zero());
215 std::fill(
T3.begin(),
T3.end(), Eigen::Matrix3d::Zero());
216 std::fill(
T2.begin(),
T2.end(), Eigen::Matrix2d::Zero());
217 std::fill(
T1.begin(),
T1.end(), Eigen::Matrix<double, 1, 1>::Zero());
228 for (
size_t i = 0; i <
T4.size(); i++)
T4[i].noalias() = A.
T4[i] * B.
T4[i];
229 for (
size_t i = 0; i <
T3.size(); i++)
T3[i].noalias() = A.
T3[i] * B.
T3[i];
230 for (
size_t i = 0; i <
T2.size(); i++)
T2[i].noalias() = A.
T2[i] * B.
T2[i];
231 for (
size_t i = 0; i <
T1.size(); i++)
T1[i].noalias() = A.
T1[i] * B.
T1[i];
243 for (
size_t i = 0; i <
T4.size(); i++)
T4[i] = A.
T4[i] * B.
T4[i];
244 for (
size_t i = 0; i <
T3.size(); i++)
T3[i] = A.
T3[i] * B.
T3[i];
245 for (
size_t i = 0; i <
T2.size(); i++)
T2[i] = A.
T2[i] * B.
T2[i];
246 for (
size_t i = 0; i <
T1.size(); i++)
T1[i] = A.
T1[i] * B.
T1[i];
289 for (
size_t i = 0; i <
T4.size(); i++)
290 T4[i].noalias() = lstm.
scale * lstm.
bas.Mat4(i);
291 for (
size_t i = 0; i <
T3.size(); i++)
292 T3[i].noalias() = lstm.
scale * lstm.
bas.Mat3(i);
293 for (
size_t i = 0; i <
T2.size(); i++)
294 T2[i].noalias() = lstm.
scale * lstm.
bas.Mat2(i);
295 for (
size_t i = 0; i <
T1.size(); i++)
296 T1[i].noalias() = lstm.
scale * lstm.
bas.Mat1(i);
325 friend std::ostream&
operator<<(std::ostream& os,
332 template <
int N>
auto&
TraMat(
size_t i)
noexcept {
333 static_assert (
N > 0 and
N < 5,
"Bad size N");
334 if constexpr (
N == 1)
return T1[i];
335 else if constexpr (
N == 2)
return T2[i];
336 else if constexpr (
N == 3)
return T3[i];
337 else if constexpr (
N == 4)
return T4[i];
341 template <
int N> [[nodiscard]]
auto&
TraMat(
size_t i)
const noexcept {
342 static_assert (
N > 0 and
N < 5,
"Bad size N");
343 if constexpr (
N == 1)
return T1[i];
344 else if constexpr (
N == 2)
return T2[i];
345 else if constexpr (
N == 3)
return T3[i];
346 else if constexpr (
N == 4)
return T4[i];
350 template <
int N> [[nodiscard]] Eigen::Matrix<Numeric, N, N>
OptDepth(
size_t i)
const noexcept {
351 return (-TraMat<N>(i).diagonal().array().log().
matrix()).asDiagonal();
375 const Eigen::Matrix<Numeric, N, 1> far,
376 const Eigen::Matrix<Numeric, N, 1> close,
377 const Eigen::Matrix<Numeric, N, N> Kfar,
378 const Eigen::Matrix<Numeric, N, N> Kclose,
381 static_assert (
N>0 and
N<5,
"Bad size N");
384 if (T(0, 0) < 0.99) {
385 const auto od = 0.5 * r * (Kfar + Kclose);
386 Eigen::Matrix<Numeric, N, 1> second = od.inverse() * ((I - (I + od) * T) * far + (od - I + T) * close);
387 if ((far[0] < close[0] and Kfar(0, 0) > Kclose(0, 0)) or (far[0] > close[0] and Kfar(0, 0) < Kclose(0, 0))) {
388 Eigen::Matrix<Numeric, N, 1> second_limit = 0.5 * (I - T) * (far + close);
393 if (second_limit[0] > second[0])
401 return 0.5 * (I - T) * (far + close);
408 const Eigen::Matrix<Numeric, N, 1> far,
409 const Eigen::Matrix<Numeric, N, 1> close,
410 const Eigen::Matrix<Numeric, N, 1>
d,
411 bool isfar)
const noexcept {
412 static_assert (
N > 0 and
N < 5,
"Bad size N");
415 const auto T = TraMat<N>(i);
416 const auto& dTdx =
dx.TraMat<
N>(i);
417 const Eigen::Matrix<Numeric, N, 1> first = 0.5 * (I - T) * (far + close);
418 if (T(0, 0) < 0.99) {
419 const auto od = OptDepth<N>(i);
420 const auto doddx =
dx.OptDepth<
N>(i);
421 const Eigen::Matrix<Numeric, N, 1> second = od.inverse() * ((I - (I + od) * T) * far + (od - I + T) * close);
422 if (first[0] > second[0]) {
424 return od.inverse() * ((I - (I + od) * T) *
d
425 - (I + od) * dTdx * far
426 + (doddx + T) * close
430 return od.inverse() * (
431 - (I + od) * dTdx * far
433 + (doddx + T) * close
437 return 0.5 * ((I - T) *
d - dTdx * (far + close));
440 return 0.5 * ((I - T) *
d - dTdx * (far + close));
469 std::vector<Eigen::Vector4d, Eigen::aligned_allocator<Eigen::Vector4d>>
R4;
470 std::vector<Eigen::Vector3d, Eigen::aligned_allocator<Eigen::Vector3d>>
R3;
471 std::vector<Eigen::Vector2d, Eigen::aligned_allocator<Eigen::Vector2d>>
R2;
472 std::vector<Eigen::Matrix<double, 1, 1>,
473 Eigen::aligned_allocator<Eigen::Matrix<double, 1, 1>>>
496 R4(std::move(rv.R4)),
497 R3(std::move(rv.R3)),
498 R2(std::move(rv.R2)),
499 R1(std::move(rv.R1)) {}
521 R4 = std::move(rv.R4);
522 R3 = std::move(rv.R3);
523 R2 = std::move(rv.R2);
524 R1 = std::move(rv.R1);
533 for (
size_t i = 0; i <
R4.size(); i++)
R4[i] = T.
Mat4(i) *
R4[i];
534 for (
size_t i = 0; i <
R3.size(); i++)
R3[i] = T.
Mat3(i) *
R3[i];
535 for (
size_t i = 0; i <
R2.size(); i++)
R2[i] = T.
Mat2(i) *
R2[i];
536 for (
size_t i = 0; i <
R1.size(); i++)
R1[i] = T.
Mat1(i) *
R1[i];
546 R4[i].noalias() = Eigen::Vector4d::Zero();
549 R3[i].noalias() = Eigen::Vector3d::Zero();
552 R2[i].noalias() = Eigen::Vector2d::Zero();
565 [[nodiscard]]
const Eigen::Vector4d&
Vec4(
size_t i)
const {
return R4[i]; }
572 [[nodiscard]]
const Eigen::Vector3d&
Vec3(
size_t i)
const {
return R3[i]; }
579 [[nodiscard]]
const Eigen::Vector2d&
Vec2(
size_t i)
const {
return R2[i]; }
586 [[nodiscard]]
const Eigen::Matrix<double, 1, 1>&
Vec1(
size_t i)
const {
return R1[i]; }
593 [[nodiscard]] Eigen::VectorXd
Vec(
size_t i)
const {
612 Eigen::Vector4d&
Vec4(
size_t i) {
return R4[i]; }
619 Eigen::Vector3d&
Vec3(
size_t i) {
return R3[i]; }
626 Eigen::Vector2d&
Vec2(
size_t i) {
return R2[i]; }
633 Eigen::Matrix<double, 1, 1>&
Vec1(
size_t i) {
return R1[i]; }
641 for (
size_t i = 0; i <
R4.size(); i++)
642 R4[i].noalias() -= 0.5 * (O1.
R4[i] + O2.
R4[i]);
643 for (
size_t i = 0; i <
R3.size(); i++)
644 R3[i].noalias() -= 0.5 * (O1.
R3[i] + O2.
R3[i]);
645 for (
size_t i = 0; i <
R2.size(); i++)
646 R2[i].noalias() -= 0.5 * (O1.
R2[i] + O2.
R2[i]);
647 for (
size_t i = 0; i <
R1.size(); i++)
648 R1[i].noalias() -= 0.5 * (O1.
R1[i] + O2.
R1[i]);
658 for (
size_t i = 0; i <
R4.size(); i++)
659 R4[i].noalias() += 0.5 * (O1.
R4[i] + O2.
R4[i]);
660 for (
size_t i = 0; i <
R3.size(); i++)
661 R3[i].noalias() += 0.5 * (O1.
R3[i] + O2.
R3[i]);
662 for (
size_t i = 0; i <
R2.size(); i++)
663 R2[i].noalias() += 0.5 * (O1.
R2[i] + O2.
R2[i]);
664 for (
size_t i = 0; i <
R1.size(); i++)
665 R1[i].noalias() += 0.5 * (O1.
R1[i] + O2.
R1[i]);
676 for (
size_t i = 0; i <
R4.size(); i++) {
679 for (
size_t i = 0; i <
R3.size(); i++) {
682 for (
size_t i = 0; i <
R2.size(); i++) {
685 for (
size_t i = 0; i <
R1.size(); i++) {
703 for (
size_t i = 0; i <
R4.size(); i++)
704 R4[i].noalias() += PiT.
Mat4(i) * (dT.
Mat4(i) * ImJ.
R4[i] + dJ.
R4[i] -
706 for (
size_t i = 0; i <
R3.size(); i++)
707 R3[i].noalias() += PiT.
Mat3(i) * (dT.
Mat3(i) * ImJ.
R3[i] + dJ.
R3[i] -
709 for (
size_t i = 0; i <
R2.size(); i++)
710 R2[i].noalias() += PiT.
Mat2(i) * (dT.
Mat2(i) * ImJ.
R2[i] + dJ.
R2[i] -
712 for (
size_t i = 0; i <
R1.size(); i++)
713 R1[i].noalias() += PiT.
Mat1(i) * (dT.
Mat1(i) * ImJ.
R1[i] + dJ.
R1[i] -
733 for (
size_t i = 0; i <
R4.size(); i++)
735 for (
size_t i = 0; i <
R3.size(); i++)
737 for (
size_t i = 0; i <
R2.size(); i++)
739 for (
size_t i = 0; i <
R1.size(); i++)
753 for (
size_t i = 0; i <
R4.size(); i++)
755 for (
size_t i = 0; i <
R3.size(); i++)
757 for (
size_t i = 0; i <
R2.size(); i++)
759 for (
size_t i = 0; i <
R1.size(); i++)
775 for (
size_t i = 0; i <
R4.size(); i++)
776 R4[i].noalias() += A.
Mat4(i) * x.
R4[i];
777 for (
size_t i = 0; i <
R3.size(); i++)
778 R3[i].noalias() += A.
Mat3(i) * x.
R3[i];
779 for (
size_t i = 0; i <
R2.size(); i++)
780 R2[i].noalias() += A.
Mat2(i) * x.
R2[i];
781 for (
size_t i = 0; i <
R1.size(); i++)
782 R1[i].noalias() += A.
Mat1(i) * x.
R1[i];
796 for (
size_t i = 0; i <
R4.size(); i++)
798 for (
size_t i = 0; i <
R3.size(); i++)
800 for (
size_t i = 0; i <
R2.size(); i++)
802 for (
size_t i = 0; i <
R1.size(); i++)
803 R4[i][0] = PiT.
Mat1(i)[0] *
818 for (
size_t i = 0; i <
R4.size(); i++)
820 for (
size_t i = 0; i <
R3.size(); i++)
822 for (
size_t i = 0; i <
R2.size(); i++)
824 for (
size_t i = 0; i <
R1.size(); i++)
839 for (
size_t i = 0; i <
R4.size(); i++)
841 for (
size_t i = 0; i <
R3.size(); i++)
843 for (
size_t i = 0; i <
R2.size(); i++)
845 for (
size_t i = 0; i <
R1.size(); i++)
856 for (
size_t i = 0; i <
R4.size(); i++) {
862 for (
size_t i = 0; i <
R3.size(); i++) {
867 for (
size_t i = 0; i <
R2.size(); i++) {
871 for (
size_t i = 0; i <
R1.size(); i++) {
902 for (
size_t i = 0; i <
R4.size(); i++)
903 for (
size_t j = 0; j < 4; j++)
M(i, j) =
R4[i](j);
904 for (
size_t i = 0; i <
R3.size(); i++)
905 for (
size_t j = 0; j < 3; j++)
M(i, j) =
R3[i](j);
906 for (
size_t i = 0; i <
R2.size(); i++)
907 for (
size_t j = 0; j < 2; j++)
M(i, j) =
R2[i](j);
908 for (
size_t i = 0; i <
R1.size(); i++)
M(i, 0) =
R1[i](0);
931 Eigen::Vector4d(
a.Kjj()[i],
a.K12()[i],
a.K13()[i],
a.K14()[i]) *
933 Eigen::Vector4d(S.
Kjj()[i], S.
K12()[i], S.
K13()[i], S.
K14()[i]);
936 Eigen::Vector4d(
a.Kjj()[i],
a.K12()[i],
a.K13()[i],
a.K14()[i]) *
942 Eigen::Vector3d(
a.Kjj()[i],
a.K12()[i],
a.K13()[i]) * B[i] +
943 Eigen::Vector3d(S.
Kjj()[i], S.
K12()[i], S.
K13()[i]);
946 Eigen::Vector3d(
a.Kjj()[i],
a.K12()[i],
a.K13()[i]) * B[i];
950 R2[i].noalias() = Eigen::Vector2d(
a.Kjj()[i],
a.K12()[i]) * B[i] +
951 Eigen::Vector2d(S.
Kjj()[i], S.
K12()[i]);
953 R2[i].noalias() = Eigen::Vector2d(
a.Kjj()[i],
a.K12()[i]) * B[i];
957 R1[i][0] =
a.Kjj()[i] * B[i] + S.
Kjj()[i];
959 R1[i][0] =
a.Kjj()[i] * B[i];
1100 const bool& jacobian_do);
1126 const Index temp_deriv_pos);
This can be used to make arrays out of anything.
A constant view of a Matrix.
A constant view of a Tensor5.
A constant view of a Vector.
Class to help with hidden temporary variables for operations of type Numeric times Class.
const base & bas
A const reference to a value.
const Numeric & scale
A const reference to a Numeric.
VectorView K13(const Index iz=0, const Index ia=0)
Vector view to K(0, 2) elements.
VectorView Kjj(const Index iz=0, const Index ia=0)
Vector view to diagonal elements.
Index NumberOfAzimuthAngles() const
The number of azimuth angles of the propagation matrix.
VectorView K12(const Index iz=0, const Index ia=0)
Vector view to K(0, 1) elements.
VectorView K14(const Index iz=0, const Index ia=0)
Vector view to K(0, 3) elements.
bool IsEmpty() const
Asks if the class is empty.
Index NumberOfZenithAngles() const
The number of zenith angles of the propagation matrix.
Stokes vector is as Propagation matrix but only has 4 possible values.
#define ARTS_ASSERT(condition,...)
Routines for setting up the jacobian.
void transform(VectorView y, double(&my_func)(double), ConstVectorView x)
A generic transform function for vectors, which can be used to implement mathematical functions opera...
NUMERIC Numeric
The type to use for all floating point numbers.
INDEX Index
The type to use for all integer numbers and indices.
Matrix matrix(const std::array< Numeric, 196 > data)
Create the square matrix from the static data.
invlib::MatrixIdentity< Matrix > Identity
invlib::Matrix< ArtsMatrix > Matrix
invlib wrapper type for ARTS matrices.
Stuff related to the propagation matrix.
Radiation Vector for Stokes dimension 1-4.
std::vector< Eigen::Matrix< double, 1, 1 >, Eigen::aligned_allocator< Eigen::Matrix< double, 1, 1 > > > R1
RadiationVector(RadiationVector &&rv) noexcept
Construct a new Radiation Vector object.
RadiationVector(Index nf=0, Index stokes=1)
Construct a new Radiation Vector object.
Eigen::Vector2d & Vec2(size_t i)
Return Vector at position.
Index Frequencies() const
Get frequency count.
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.
Eigen::Matrix< double, 1, 1 > & Vec1(size_t i)
Return Vector at position.
void SetZero(size_t i)
Set Radiation Vector to Zero at position.
void rem_avg(const RadiationVector &O1, const RadiationVector &O2)
Remove the average of two other RadiationVector from *this.
std::vector< Eigen::Vector2d, Eigen::aligned_allocator< Eigen::Vector2d > > R2
void add_avg(const RadiationVector &O1, const RadiationVector &O2)
Add the average of two other RadiationVector to *this.
RadiationVector & operator=(RadiationVector &&rv) noexcept
Assign old radiation vector to this.
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.
std::vector< Eigen::Vector4d, Eigen::aligned_allocator< Eigen::Vector4d > > R4
const Eigen::Vector2d & Vec2(size_t i) const
Return Vector at position.
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.
Eigen::Vector4d & Vec4(size_t i)
Return Vector at position.
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 ConstMatrixView &M)
Set *this from matrix.
RadiationVector & operator=(const RadiationVector &rv)=default
Assign old radiation vector to this.
std::vector< Eigen::Vector3d, Eigen::aligned_allocator< Eigen::Vector3d > > R3
Eigen::Vector3d & Vec3(size_t i)
Return Vector at position.
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.
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=(const LazyScale< TransmissionMatrix > &lstm)
Assign lazily.
TransmissionMatrix(TransmissionMatrix &&tm) noexcept
Construct a new Transmission Matrix object.
Eigen::Matrix3d & Mat3(size_t i)
Get Matrix at position.
std::vector< Eigen::Matrix3d, Eigen::aligned_allocator< Eigen::Matrix3d > > T3
void mul_aliased(const TransmissionMatrix &A, const TransmissionMatrix &B)
Set this to a multiple of A by B.
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.
Eigen::Matrix4d & Mat4(size_t i)
Get Matrix at position.
const Eigen::Matrix2d & Mat2(size_t i) const
Get Matrix at position.
TransmissionMatrix(const TransmissionMatrix &tm)=default
Construct a new Transmission Matrix object.
std::vector< Eigen::Matrix< double, 1, 1 >, Eigen::aligned_allocator< Eigen::Matrix< double, 1, 1 > > > T1
TransmissionMatrix & operator=(TransmissionMatrix &&tm) noexcept
Move operator.
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.
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.
friend std::istream & operator>>(std::istream &data, TransmissionMatrix &tm)
Input operator.
const Eigen::Matrix3d & Mat3(size_t i) const
Get Matrix at position.
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.
std::vector< Eigen::Matrix2d, Eigen::aligned_allocator< Eigen::Matrix2d > > T2
TransmissionMatrix & operator=(const TransmissionMatrix &tm)=default
Assignment operator.
TransmissionMatrix & operator+=(const LazyScale< TransmissionMatrix > &lstm)
Assign to *this lazily.
void setZero()
Set to zero matrix.
Eigen::Matrix< double, 1, 1 > & Mat1(size_t i)
Get Matrix at position.
TransmissionMatrix(Index nf=0, Index stokes=1)
Construct a new Transmission Matrix object.
Eigen::Matrix2d & Mat2(size_t i)
Get Matrix at position.
std::vector< Eigen::Matrix4d, Eigen::aligned_allocator< Eigen::Matrix4d > > T4
void stepwise_source(RadiationVector &J, ArrayOfRadiationVector &dJ, 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.
RadiativeTransferSolver
Intended to hold various forward solvers.
ArrayOfTransmissionMatrix bulk_backscatter(const ConstTensor5View &Pe, const ConstMatrixView &pnd)
Bulk back-scattering.
BackscatterSolver
Intended to hold various backscatter solvers.
@ CommutativeTransmission
CumulativeTransmission
Intended to hold various ways to accumulate the transmission matrix.
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.
std::istream & operator>>(std::istream &is, TransmissionMatrix &tm)
Input operator.
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.
std::ostream & operator<<(std::ostream &os, const TransmissionMatrix &tm)
Output operator.
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.