diff --git a/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/GaussianSource.h b/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/GaussianSource.h
index 50616ae51c04bbf82db3ec41fdf12e6b1965ffab..6d4c4fe3202a52da2a919ce7e596026dce86357f 100644
--- a/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/GaussianSource.h
+++ b/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/GaussianSource.h
@@ -46,7 +46,7 @@ public:
GaussianSource(const string &name, const Expr &ra, const Expr &dec,
const Expr &I, const Expr &Q, const Expr &U, const Expr &V,
- const Expr &si, const Expr &major, const Expr &minor, const Expr &phi);
+ const Expr &major, const Expr &minor, const Expr &phi);
const Expr &getI() const
{ return itsI; }
@@ -56,8 +56,6 @@ public:
{ return itsU; }
const Expr &getV() const
{ return itsV; }
- const Expr &getSpectralIndex() const
- { return itsSpectralIndex; }
const Expr &getMajor() const
{ return itsMajor; }
const Expr &getMinor() const
@@ -73,7 +71,6 @@ private:
Expr itsMajor;
Expr itsMinor;
Expr itsPhi;
- Expr itsSpectralIndex;
};
// @}
diff --git a/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/PointSource.h b/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/PointSource.h
index 5135a9fd0a4a2ceb1582a1da36b88346b72db1e7..32d042ce8204accb1f78a3a05c64a6073fcc507c 100644
--- a/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/PointSource.h
+++ b/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/PointSource.h
@@ -43,9 +43,8 @@ public:
typedef shared_ptr<PointSource> Pointer;
typedef shared_ptr<const PointSource> ConstPointer;
- PointSource(const string& name, const Expr &ra, const Expr &dec,
- const Expr &I, const Expr &Q, const Expr &U, const Expr &V,
- const Expr &si);
+ PointSource(const string &name, const Expr &ra, const Expr &dec,
+ const Expr &I, const Expr &Q, const Expr &U, const Expr &V);
const Expr &getI() const
{ return itsI; }
@@ -55,15 +54,12 @@ public:
{ return itsU; }
const Expr &getV() const
{ return itsV; }
- const Expr &getSpectralIndex() const
- { return itsSpectralIndex; }
private:
Expr itsI;
Expr itsQ;
Expr itsU;
Expr itsV;
- Expr itsSpectralIndex;
};
// @}
diff --git a/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/SpectralIndex.h b/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/SpectralIndex.h
index 707d824aa3d8df6029934af905d84a8fdbebc922..c8ad90d5f30b1bf06654c5167a759bc02e59912c 100644
--- a/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/SpectralIndex.h
+++ b/CEP/Calibration/BBSKernel/include/BBSKernel/Expr/SpectralIndex.h
@@ -42,17 +42,19 @@ class SpectralIndex: public ExprRep
{
public:
template <typename T_Iterator>
- SpectralIndex(const Expr &refFreq, T_Iterator first, T_Iterator last);
+ SpectralIndex(const Expr &refFreq, const Expr &refStokes, T_Iterator first,
+ T_Iterator last);
virtual Matrix getResultValue(const Request &request,
const std::vector<const Matrix*> &args);
};
template <typename T_Iterator>
-SpectralIndex::SpectralIndex(const Expr &refFreq, T_Iterator first,
- T_Iterator last)
+SpectralIndex::SpectralIndex(const Expr &refFreq, const Expr &refStokes,
+ T_Iterator first, T_Iterator last)
{
addChild(refFreq);
+ addChild(refStokes);
while(first != last)
{
addChild(*first++);
diff --git a/CEP/Calibration/BBSKernel/include/BBSKernel/Model.h b/CEP/Calibration/BBSKernel/include/BBSKernel/Model.h
index 7f76c7aa2115058d7558a691ef6d917b1bfd2836..8b846468e945b8aa27f1897fdd8b66580269a313 100644
--- a/CEP/Calibration/BBSKernel/include/BBSKernel/Model.h
+++ b/CEP/Calibration/BBSKernel/include/BBSKernel/Model.h
@@ -103,7 +103,7 @@ private:
Source::Pointer makeSource(const SourceInfo &source);
- Expr makeSpectralIndex(const SourceInfo &source);
+ Expr makeSpectralIndex(const SourceInfo &source, const string &stokesParm);
void makeStationUvw();
diff --git a/CEP/Calibration/BBSKernel/src/Expr/GaussianCoherence.cc b/CEP/Calibration/BBSKernel/src/Expr/GaussianCoherence.cc
index 16eea22a0608f59c5a4519b65443f5c5f6e7db99..2ce6acf1be0a6c08da5c008d9966a8afae7139ed 100644
--- a/CEP/Calibration/BBSKernel/src/Expr/GaussianCoherence.cc
+++ b/CEP/Calibration/BBSKernel/src/Expr/GaussianCoherence.cc
@@ -53,7 +53,6 @@ GaussianCoherence::GaussianCoherence(const GaussianSource::ConstPointer &source,
addChild(source->getQ());
addChild(source->getU());
addChild(source->getV());
- addChild(source->getSpectralIndex());
addChild(source->getMajor());
addChild(source->getMinor());
addChild(source->getPhi());
@@ -63,7 +62,7 @@ JonesResult GaussianCoherence::getJResult(const Request &request)
{
enum ChildExprIndex
{
- IN_I, IN_Q, IN_U, IN_V, IN_SI, IN_MAJOR, IN_MINOR, IN_PHI
+ IN_I, IN_Q, IN_U, IN_V, IN_MAJOR, IN_MINOR, IN_PHI
};
// Evaluate source parameters.
@@ -74,7 +73,6 @@ JonesResult GaussianCoherence::getJResult(const Request &request)
const Result &qk = getChild(IN_Q).getResultSynced(request, qkBuf);
const Result &uk = getChild(IN_U).getResultSynced(request, ukBuf);
const Result &vk = getChild(IN_V).getResultSynced(request, vkBuf);
- const Result &si = getChild(IN_SI).getResultSynced(request, siBuf);
const Result &major = getChild(IN_MAJOR).getResultSynced(request, majorBuf);
const Result &minor = getChild(IN_MINOR).getResultSynced(request, minorBuf);
const Result &phi = getChild(IN_PHI).getResultSynced(request, phiBuf);
@@ -118,68 +116,61 @@ JonesResult GaussianCoherence::getJResult(const Request &request)
Matrix vBaseline = vStation2.getValue() - vStation1.getValue();
// Compute spatial coherence for a gaussian source.
- Matrix coherence = computeCoherence(request, uBaseline, vBaseline,
- major.getValue(), minor.getValue(), phi.getValue());
+ Matrix coherence_2 = computeCoherence(request, uBaseline, vBaseline,
+ major.getValue(), minor.getValue(), phi.getValue()) * 0.5;
// Compute main result.
- Matrix uv_2 = tocomplex(uk.getValue(), vk.getValue()) * 0.5;
+ Matrix uv = tocomplex(uk.getValue(), vk.getValue());
- resXX.setValue(si.getValue() * coherence * (ik.getValue() + qk.getValue())
- * 0.5);
- resXY.setValue(uv_2 * coherence);
- resYX.setValue(conj(uv_2) * coherence);
- resYY.setValue(si.getValue() * coherence * (ik.getValue() - qk.getValue())
- * 0.5);
+ resXX.setValue(coherence_2 * (ik.getValue() + qk.getValue()));
+ resXY.setValue(coherence_2 * uv);
+ resYX.setValue(coherence_2 * conj(uv));
+ resYY.setValue(coherence_2 * (ik.getValue() - qk.getValue()));
// Compute the perturbed values.
enum PValues
{
- PV_I, PV_Q, PV_U, PV_V, PV_SI, PV_MAJOR, PV_MINOR, PV_PHI, N_PValues
+ PV_I, PV_Q, PV_U, PV_V, PV_MAJOR, PV_MINOR, PV_PHI, N_PValues
};
- const Result *pvSet[N_PValues] = {&ik, &qk, &uk, &vk, &si, &major, &minor,
- &phi};
+ const Result *pvSet[N_PValues] = {&ik, &qk, &uk, &vk, &major, &minor, &phi};
PValueSetIterator<N_PValues> pvIter(pvSet);
while(!pvIter.atEnd())
{
- bool evalIQ = pvIter.hasPValue(PV_I) || pvIter.hasPValue(PV_Q)
- || pvIter.hasPValue(PV_SI);
+ bool evalIQ = pvIter.hasPValue(PV_I) || pvIter.hasPValue(PV_Q);
bool evalUV = pvIter.hasPValue(PV_U) || pvIter.hasPValue(PV_V);
bool evalCoh = pvIter.hasPValue(PV_MAJOR) || pvIter.hasPValue(PV_MINOR)
|| pvIter.hasPValue(PV_PHI);
- Matrix pvUV_2(uv_2);
- Matrix pvCoh(coherence);
+ Matrix pvUV(uv);
+ Matrix pvCoh_2(coherence_2);
if(evalUV)
{
- pvUV_2 = tocomplex(pvIter.value(PV_U), pvIter.value(PV_V)) * 0.5;
+ pvUV = tocomplex(pvIter.value(PV_U), pvIter.value(PV_V));
}
if(evalCoh)
{
- pvCoh = computeCoherence(request, uBaseline, vBaseline,
+ pvCoh_2 = computeCoherence(request, uBaseline, vBaseline,
pvIter.value(PV_MAJOR), pvIter.value(PV_MINOR),
- pvIter.value(PV_PHI));
+ pvIter.value(PV_PHI)) * 0.5;
}
if(evalIQ || evalCoh)
{
const Matrix &pvI = pvIter.value(PV_I);
const Matrix &pvQ = pvIter.value(PV_Q);
- const Matrix &pvSI = pvIter.value(PV_SI);
- resXX.setPerturbedValue(pvIter.key(), pvSI * pvCoh * (pvI + pvQ)
- * 0.5);
- resYY.setPerturbedValue(pvIter.key(), pvSI * pvCoh * (pvI - pvQ)
- * 0.5);
+ resXX.setPerturbedValue(pvIter.key(), pvCoh_2 * (pvI + pvQ));
+ resYY.setPerturbedValue(pvIter.key(), pvCoh_2 * (pvI - pvQ));
}
if(evalUV || evalCoh)
{
- resXY.setPerturbedValue(pvIter.key(), pvUV_2 * pvCoh);
- resYX.setPerturbedValue(pvIter.key(), conj(pvUV_2) * pvCoh);
+ resXY.setPerturbedValue(pvIter.key(), pvCoh_2 * pvUV);
+ resYX.setPerturbedValue(pvIter.key(), pvCoh_2 * conj(pvUV));
}
pvIter.next();
diff --git a/CEP/Calibration/BBSKernel/src/Expr/GaussianSource.cc b/CEP/Calibration/BBSKernel/src/Expr/GaussianSource.cc
index 9fc290ccae39e79fe077a5b338ba280dd40a5a8f..291b85d93961adc16db432276c421bdf878243bf 100644
--- a/CEP/Calibration/BBSKernel/src/Expr/GaussianSource.cc
+++ b/CEP/Calibration/BBSKernel/src/Expr/GaussianSource.cc
@@ -30,14 +30,13 @@ namespace BBS
GaussianSource::GaussianSource(const string &name,
const Expr &ra, const Expr &dec,
- const Expr &I, const Expr &Q, const Expr &U, const Expr &V, const Expr &si,
+ const Expr &I, const Expr &Q, const Expr &U, const Expr &V,
const Expr &major, const Expr &minor, const Expr &phi)
: Source(name, ra, dec),
itsI(I),
itsQ(Q),
itsU(U),
itsV(V),
- itsSpectralIndex(si),
itsMajor(major),
itsMinor(minor),
itsPhi(phi)
diff --git a/CEP/Calibration/BBSKernel/src/Expr/PointCoherence.cc b/CEP/Calibration/BBSKernel/src/Expr/PointCoherence.cc
index 322d5f0a0c82018764dc00652c3549b0e410128d..634649348c80b4061bb4c428829acca8cc14e16b 100644
--- a/CEP/Calibration/BBSKernel/src/Expr/PointCoherence.cc
+++ b/CEP/Calibration/BBSKernel/src/Expr/PointCoherence.cc
@@ -40,7 +40,6 @@ PointCoherence::PointCoherence(const PointSource::ConstPointer &source)
addChild(source->getQ());
addChild(source->getU());
addChild(source->getV());
- addChild(source->getSpectralIndex());
}
@@ -61,39 +60,36 @@ JonesResult PointCoherence::getJResult(const Request &request)
const Result& qk = getChild(1).getResultSynced(request, qkBuf);
const Result& uk = getChild(2).getResultSynced(request, ukBuf);
const Result& vk = getChild(3).getResultSynced(request, vkBuf);
- const Result& si = getChild(4).getResultSynced(request, siBuf);
// Compute main value.
- Matrix uvk_2 = tocomplex(uk.getValue(), vk.getValue()) * 0.5;
- resXX.setValue(si.getValue() * (ik.getValue() + qk.getValue()) * 0.5);
+ Matrix uvk_2 = 0.5 * tocomplex(uk.getValue(), vk.getValue());
+ resXX.setValue(0.5 * (ik.getValue() + qk.getValue()));
resXY.setValue(uvk_2);
resYX.setValue(conj(uvk_2));
- resYY.setValue(si.getValue() * (ik.getValue() - qk.getValue()) * 0.5);
+ resYY.setValue(0.5 * (ik.getValue() - qk.getValue()));
// Compute perturbed values.
enum PValues
- { PV_I, PV_Q, PV_U, PV_V, PV_SI, N_PValues };
+ { PV_I, PV_Q, PV_U, PV_V, N_PValues };
- const Result *pvSet[N_PValues] = {&ik, &qk, &uk, &vk, &si};
+ const Result *pvSet[N_PValues] = {&ik, &qk, &uk, &vk};
PValueSetIterator<N_PValues> pvIter(pvSet);
while(!pvIter.atEnd())
{
- if(pvIter.hasPValue(PV_I) || pvIter.hasPValue(PV_Q)
- || pvIter.hasPValue(PV_SI))
+ if(pvIter.hasPValue(PV_I) || pvIter.hasPValue(PV_Q))
{
const Matrix &pvI = pvIter.value(PV_I);
const Matrix &pvQ = pvIter.value(PV_Q);
- const Matrix &pvSI = pvIter.value(PV_SI);
- resXX.setPerturbedValue(pvIter.key(), pvSI * (pvI + pvQ) * 0.5);
- resYY.setPerturbedValue(pvIter.key(), pvSI * (pvI - pvQ) * 0.5);
+ resXX.setPerturbedValue(pvIter.key(), 0.5 * (pvI + pvQ));
+ resYY.setPerturbedValue(pvIter.key(), 0.5 * (pvI - pvQ));
}
if(pvIter.hasPValue(PV_U) || pvIter.hasPValue(PV_V))
{
const Matrix &pvU = pvIter.value(PV_U);
const Matrix &pvV = pvIter.value(PV_V);
- Matrix uvkp_2 = tocomplex(pvU, pvV) * 0.5;
+ Matrix uvkp_2 = 0.5 * tocomplex(pvU, pvV);
resXY.setPerturbedValue(pvIter.key(), uvkp_2);
resYX.setPerturbedValue(pvIter.key(), conj(uvkp_2));
}
diff --git a/CEP/Calibration/BBSKernel/src/Expr/PointSource.cc b/CEP/Calibration/BBSKernel/src/Expr/PointSource.cc
index 47047dfe39502c1d58025bee4ed83065df46ee83..f16f4eae11f507cbf7a4455c1cc918066e135a7d 100644
--- a/CEP/Calibration/BBSKernel/src/Expr/PointSource.cc
+++ b/CEP/Calibration/BBSKernel/src/Expr/PointSource.cc
@@ -28,14 +28,13 @@ namespace LOFAR
namespace BBS
{
-PointSource::PointSource(const string& name, const Expr &ra, const Expr &dec,
- const Expr &I, const Expr &Q, const Expr &U, const Expr &V, const Expr &si)
+PointSource::PointSource(const string &name, const Expr &ra, const Expr &dec,
+ const Expr &I, const Expr &Q, const Expr &U, const Expr &V)
: Source(name, ra, dec),
itsI(I),
itsQ(Q),
itsU(U),
- itsV(V),
- itsSpectralIndex(si)
+ itsV(V)
{
}
diff --git a/CEP/Calibration/BBSKernel/src/Expr/SpectralIndex.cc b/CEP/Calibration/BBSKernel/src/Expr/SpectralIndex.cc
index 66597f9dd2b521173ba84226923f9b9699ac59e9..e33f7324f958824c339c7b2253ae16c845e537f9 100644
--- a/CEP/Calibration/BBSKernel/src/Expr/SpectralIndex.cc
+++ b/CEP/Calibration/BBSKernel/src/Expr/SpectralIndex.cc
@@ -36,11 +36,12 @@ Matrix SpectralIndex::getResultValue(const Request &request,
ASSERT(args.size() == getChildCount());
// Special case: No coefficients or a single real coefficient with value
- // zero.
- if(args.size() == 1 || (args.size() == 2 && !args.back()->isArray()
+ // zero. Return the Stokes parameter at the reference frequency in this
+ // case.
+ if(args.size() == 2 || (args.size() == 3 && !args.back()->isArray()
&& !args.back()->isComplex() && args.back()->getDouble() == 0.0))
{
- return Matrix(1.0);
+ return *args[1];
}
// Create a 2-D Matrix that contains the frequency for each sample. (We
@@ -59,8 +60,8 @@ Matrix SpectralIndex::getResultValue(const Request &request,
}
}
- // Compute flux scale factor as:
- // (v / v0) ^ (-1.0 * [c0 + c1 * log(v / v0) + c2 * log(v / v0)^2 + ...])
+ // Compute scale factor as:
+ // (v / v0) ^ (c0 + c1 * log(v / v0) + c2 * log(v / v0)^2 + ...)
// Where v is the frequency and v0 is the reference frequency.
// Compute log(v / v0).
@@ -69,13 +70,14 @@ Matrix SpectralIndex::getResultValue(const Request &request,
// Compute c0 + log(v / v0) * c1 + log(v / v0)^2 * c2 + ... using Horner's
// rule.
Matrix exponent = *args.back();
- for(unsigned int i = args.size() - 2; i >= 1; --i)
+ for(unsigned int i = args.size() - 2; i >= 2; --i)
{
exponent = exponent * base + *args[i];
}
- // Compute (v / v0) ^ -exponent.
- return exp(base * (-exponent));
+ // Compute I0 * (v / v0) ^ exponent, where I0 is the value of the Stokes
+ // parameter at the reference frequency.
+ return *args[1] * exp(base * exponent);
}
} //# namespace BBS
diff --git a/CEP/Calibration/BBSKernel/src/Model.cc b/CEP/Calibration/BBSKernel/src/Model.cc
index fbe5a68c84699bafc76acc5df8310473d465eacf..2e0b8c9c1d4b34acc865c96f6f43fd24bfba3145 100644
--- a/CEP/Calibration/BBSKernel/src/Model.cc
+++ b/CEP/Calibration/BBSKernel/src/Model.cc
@@ -726,14 +726,13 @@ Source::Pointer Model::makeSource(const SourceInfo &source)
{
Expr ra = makeExprParm(SKY, "Ra:" + source.getName());
Expr dec = makeExprParm(SKY, "Dec:" + source.getName());
- Expr i = makeExprParm(SKY, "I:" + source.getName());
+ Expr i = makeSpectralIndex(source, "I");
Expr q = makeExprParm(SKY, "Q:" + source.getName());
Expr u = makeExprParm(SKY, "U:" + source.getName());
Expr v = makeExprParm(SKY, "V:" + source.getName());
- Expr spectral = makeSpectralIndex(source);
return PointSource::Pointer(new PointSource(source.getName(),
- ra, dec, i, q, u, v, spectral));
+ ra, dec, i, q, u, v));
}
break;
@@ -741,17 +740,16 @@ Source::Pointer Model::makeSource(const SourceInfo &source)
{
Expr ra = makeExprParm(SKY, "Ra:" + source.getName());
Expr dec = makeExprParm(SKY, "Dec:" + source.getName());
- Expr i = makeExprParm(SKY, "I:" + source.getName());
+ Expr i = makeSpectralIndex(source, "I");
Expr q = makeExprParm(SKY, "Q:" + source.getName());
Expr u = makeExprParm(SKY, "U:" + source.getName());
Expr v = makeExprParm(SKY, "V:" + source.getName());
Expr maj = makeExprParm(SKY, "Major:" + source.getName());
Expr min = makeExprParm(SKY, "Minor:" + source.getName());
Expr phi = makeExprParm(SKY, "Phi:" + source.getName());
- Expr spectral = makeSpectralIndex(source);
return GaussianSource::Pointer(new GaussianSource(source.getName(),
- ra, dec, i, q, u, v, spectral, maj, min, phi));
+ ra, dec, i, q, u, v, maj, min, phi));
}
break;
@@ -763,14 +761,17 @@ Source::Pointer Model::makeSource(const SourceInfo &source)
return Source::Pointer();
}
-Expr Model::makeSpectralIndex(const SourceInfo &source)
+Expr Model::makeSpectralIndex(const SourceInfo &source,
+ const string &stokesParm)
{
unsigned int degree =
static_cast<unsigned int>(ParmManager::instance().getDefaultValue(SKY,
"SpectralIndexDegree:" + source.getName()));
- Expr reference = makeExprParm(SKY, "ReferenceFrequency:"
- + source.getName());
+ // Reference frequency.
+ Expr refFreq = makeExprParm(SKY, "ReferenceFrequency:" + source.getName());
+ // Stokes parameter value at the reference frequency.
+ Expr refStokes = makeExprParm(SKY, stokesParm + ":" + source.getName());
vector<Expr> coefficients;
coefficients.reserve(degree + 1);
@@ -781,7 +782,7 @@ Expr Model::makeSpectralIndex(const SourceInfo &source)
coefficients.push_back(makeExprParm(SKY, name.str()));
}
- return Expr(new SpectralIndex(reference, coefficients.begin(),
+ return Expr(new SpectralIndex(refFreq, refStokes, coefficients.begin(),
coefficients.end()));
}