26 #ifndef GRALE_BACKPROJECTMATRIXNEW_H
28 #define GRALE_BACKPROJECTMATRIXNEW_H
34 #include "graleconfig.h"
35 #include "deflectionmatrix.h"
41 class ImagesDataExtended;
42 class GravitationalLens;
52 bool startInit(
double z_d,
double D_d, DeflectionMatrix *pDeflectionMatrix,
53 const std::vector<ImagesDataExtended *> &images,
54 const std::vector<bool> &useDeflections,
55 const std::vector<bool> &useDerivatives,
56 const std::vector<bool> &usePotentials,
58 bool storeOriginalIntensities,
59 bool storeOriginalTimeDelays,
60 bool storeOriginalShearInfo,
64 void storeDeflectionMatrixResults();
65 void calculate(
float scaleFactor,
float massSheetFactor = 0);
66 void calculateInverseMagnifications() { calculateInverseMagnifications(m_trueFlags); }
67 void calculateShearComponents() { calculateShearComponents(m_trueFlags); }
68 void calculateConvergence() { calculateConvergence(m_trueFlags); }
69 void calculateInverseMagnifications(
const std::vector<bool> &sourceMask);
70 void calculateShearComponents(
const std::vector<bool> &sourceMask);
71 void calculateConvergence(
const std::vector<bool> &sourceMask);
73 double getAngularScale()
const {
return m_angularScale; }
74 double getMassSheetScale()
const {
return m_massSheetScale; }
76 const Vector2D<float> *getBetas(
int sourceNumber)
const {
return &(m_betas[sourceNumber][0]); }
77 const Vector2D<float> *getBetas(
int sourceNumber,
int imageNumber)
const {
return &(m_betas[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
78 const Vector2D<float> *getThetas(
int sourceNumber)
const {
return &(m_thetas[sourceNumber][0]); }
79 const Vector2D<float> *getThetas(
int sourceNumber,
int imageNumber)
const {
return &(m_thetas[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
80 const float *getDerivativesXX(
int sourceNumber)
const {
return &(m_axx[sourceNumber][0]); }
81 const float *getDerivativesXX(
int sourceNumber,
int imageNumber)
const {
return &(m_axx[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
82 const float *getDerivativesYY(
int sourceNumber)
const {
return &(m_ayy[sourceNumber][0]); }
83 const float *getDerivativesYY(
int sourceNumber,
int imageNumber)
const {
return &(m_ayy[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
84 const float *getDerivativesXY(
int sourceNumber)
const {
return &(m_axy[sourceNumber][0]); }
85 const float *getDerivativesXY(
int sourceNumber,
int imageNumber)
const {
return &(m_axy[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
86 const float *getInverseMagnifications(
int sourceNumber)
const {
return &(m_inverseMagnifications[sourceNumber][0]); }
87 const float *getInverseMagnifications(
int sourceNumber,
int imageNumber)
const {
return &(m_inverseMagnifications[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
88 const float *getShearComponents1(
int sourceNumber)
const {
return &(m_shearComponent1[sourceNumber][0]); }
89 const float *getShearComponents1(
int sourceNumber,
int imageNumber)
const {
return &(m_shearComponent1[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
90 const float *getShearComponents2(
int sourceNumber)
const {
return &(m_axy[sourceNumber][0]); }
91 const float *getShearComponents2(
int sourceNumber,
int imageNumber)
const {
return &(m_axy[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
92 const float *getConvergence(
int sourceNumber)
const {
return &(m_convergence[sourceNumber][0]); }
93 const float *getConvergence(
int sourceNumber,
int imageNumber)
const {
return &(m_convergence[sourceNumber][m_offsets[sourceNumber][imageNumber]]); }
94 float getTimeDelay(
int sourceNumber,
int imageNumber,
int pointNumber,
Vector2D<float> beta)
const;
96 DeflectionMatrix *m_pDeflectionMatrix;
97 bool m_init, m_initializing;
99 std::vector<std::vector<int> > m_deflectionIndices;
100 std::vector<std::vector<int> > m_derivativeIndices;
101 std::vector<std::vector<int> > m_potentialIndices;
103 std::vector<std::vector<Vector2D<double> > > m_originalPoints;
106 std::vector<std::vector<float> > m_baseAxx, m_baseAyy, m_baseAxy;
107 std::vector<std::vector<float> > m_basePotentials;
109 std::vector<std::vector<Vector2D<float> > > m_baseAlphas;
110 std::vector<std::vector<Vector2D<double> > > m_baseAlphasUnscaled;
112 std::vector<std::vector<Vector2D<float> > > m_baseBetas;
113 std::vector<std::vector<Vector2D<double> > > m_baseBetasUnscaled;
114 #endif // AUTODISTFRAC
115 std::vector<std::vector<double> > m_basePotentialsUnscaled;
118 double m_angularScale;
119 double m_massSheetScale;
120 float m_timeDelayScale;
122 std::vector<std::vector<Vector2D<float> > > m_thetas;
123 std::vector<std::vector<Vector2D<float> > > m_subDeflectionAngles;
124 std::vector<std::vector<float> > m_subDeflectionDerivatives[3];
125 std::vector<std::vector<float> > m_subPotentialValues;
127 std::vector<std::vector<Vector2D<float> > > m_betas;
128 std::vector<std::vector<float> > m_axx, m_ayy, m_axy;
129 std::vector<std::vector<float> > m_potentials;
130 std::vector<std::vector<float> > m_inverseMagnifications;
131 std::vector<std::vector<float> > m_shearComponent1;
132 std::vector<std::vector<float> > m_convergence;
134 std::vector<std::vector<Vector2D<float> > > m_sheetThetas;
135 std::vector<std::vector<float> > m_sheetPotentials;
138 std::vector<float> m_tmpBuffer;
139 std::vector<float> m_oneVector;
140 std::vector<bool> m_trueFlags;
143 inline float BackProjectMatrixNew::getTimeDelay(
int sourceNumber,
int imageNumber,
int pointNumber,
Vector2D<float> beta)
const
145 Vector2D<float> diff = beta - m_thetas[sourceNumber][m_offsets[sourceNumber][imageNumber]+pointNumber];
147 return m_timeDelayScale*(0.5f*diff.getLengthSquared() - m_potentials[sourceNumber][m_offsets[sourceNumber][imageNumber]+pointNumber])/m_distanceFractions[sourceNumber];
152 #endif // GRALE_BACKPROJECTMATRIXNEW_H