diff --git a/include/Directions.h b/include/Directions.h
index f172324e00e452c3b90e0cdab430d72a716e5831..58a982be7f798bc5009799da02af4ef6a77fc5e9 100644
--- a/include/Directions.h
+++ b/include/Directions.h
@@ -51,6 +51,8 @@ public:
     dec.clear();
   }
 
+  Direction operator[](size_t i) const { return Direction(ra[i], dec[i]); }
+
   size_t Size() const { return ra.size(); }
 
   SmartVector<double> ra;
diff --git a/include/GaussianSource.h b/include/GaussianSource.h
index e12c799c8ff677429a171bad82015c0f25406544..17065a5d43c0276cb9bebc0fd09683d9c1e9f2ac 100644
--- a/include/GaussianSource.h
+++ b/include/GaussianSource.h
@@ -6,7 +6,7 @@
 #ifndef DPPP_GAUSSIANSOURCE_H
 #define DPPP_GAUSSIANSOURCE_H
 
-#include "PointSource.h"
+#include <PointSource.h>
 
 namespace dp3 {
 namespace base {
@@ -15,7 +15,7 @@ namespace base {
 
 /// @{
 
-class GaussianSource : public PointSource {
+class GaussianSource : public dp3::base::PointSource {
 public:
   typedef std::shared_ptr<GaussianSource> Ptr;
   typedef std::shared_ptr<const GaussianSource> ConstPtr;
@@ -23,6 +23,9 @@ public:
   GaussianSource(const Direction &direction);
   GaussianSource(const Direction &direction, const Stokes &stokes,
                  size_t beam_id = 0);
+  GaussianSource(const Direction &direction, const Spectrum &spectrum,
+                 double position_angle, bool is_position_angle_absolute,
+                 double minor_axis, double major_axis, size_t beam_id = 0);
 
   /// Set position angle in radians. The position angle is the smallest angle
   /// between the major axis and North, measured positively North over East.
diff --git a/include/GaussianSourceCollection.h b/include/GaussianSourceCollection.h
index 9c3e2b53498ef25c9fcec214f060b67b2dc62dc7..3c0c26f98062c837de665732657bc75f022bef85 100644
--- a/include/GaussianSourceCollection.h
+++ b/include/GaussianSourceCollection.h
@@ -38,6 +38,23 @@ public:
     PointSourceCollection::Reserve(new_size);
   }
 
+  GaussianSource operator[](size_t i) const {
+    return GaussianSource(direction_vector[i], spectrums[i], beam_id[i],
+                          position_angle[i], major_axis[i], minor_axis[i],
+                          position_angle_is_absolute[i]);
+  }
+
+  std::unique_ptr<GaussianSourceCollection> SelectBeamID(size_t beam_id) {
+    auto selected = std::make_unique<GaussianSourceCollection>();
+
+    for (size_t i = 0; i < Size(); ++i) {
+      if (this->beam_id[i] == beam_id) {
+        selected->Add(operator[](i));
+      }
+    }
+    return std::move(selected);
+  }
+
   SmartVector<double> position_angle;
   SmartVector<double> major_axis;
   SmartVector<double> minor_axis;
diff --git a/include/PointSource.h b/include/PointSource.h
index 154503f7089d87204265b4b7ae89eb8e3433a60a..0b41548e0cbf26a2b12cd60df2486e676bbcc26d 100644
--- a/include/PointSource.h
+++ b/include/PointSource.h
@@ -10,7 +10,6 @@
 #include "Direction.h"
 #include "Spectrum.h"
 #include "Stokes.h"
-#include "StokesVector.h"
 #include <memory>
 #include <set>
 namespace dp3 {
@@ -28,6 +27,8 @@ public:
 
   PointSource(const Direction &direction, const Stokes &stokes,
               const size_t beam_id = 0);
+  PointSource(const Direction &direction, const Spectrum &stokes,
+              const size_t beam_id = 0);
 
   const Direction &GetDirection() const { return direction_; }
   void SetDirection(const Direction &position);
@@ -37,8 +38,7 @@ public:
   const Spectrum &GetSpectrum() const { return spectrum_; }
 
   Stokes GetStokes(double freq) const;
-  const Stokes &GetStokes() const { return stokes_; }
-  const StokesVector &GetStokesVector() const { return stokes_vector_; }
+  const Stokes &GetStokes() const { return spectrum_.GetReferenceFlux(); }
 
   const size_t &GetBeamId() const { return beam_id_; }
   void SetBeamId(size_t beam_id) { beam_id_ = beam_id; }
@@ -54,19 +54,7 @@ public:
 protected:
   Direction direction_;
   Spectrum spectrum_;
-  StokesVector stokes_vector_;
-  Stokes stokes_;
   size_t beam_id_;
-  // FIXME: remove the following declerations.
-  // They are not needed anymore, but are kept for
-  // backward compatibility for the tests.
-  double reference_frequency_;
-  std::vector<double> spectral_terms_;
-  double polarizated_fraction_;
-  double polarization_angle_;
-  double rotation_measure_;
-  bool has_rotation_measure_;
-  bool has_logarithmic_si_;
 };
 
 /// @}
@@ -77,18 +65,10 @@ void PointSource::SetSpectralTerms(double refFreq, bool isLogarithmic, T first,
   // FIXME: the following four statements can be removed later.
   // They are not needed anymore, but are kept for
   // backward compatibility for the tests.
-  reference_frequency_ = refFreq;
-  has_logarithmic_si_ = isLogarithmic;
-
-  spectral_terms_.clear();
-  spectral_terms_.insert(spectral_terms_.begin(), first, last);
-  auto new_terms = xt::adapt(std::vector<double>(first, last));
-
-  spectrum_.SetSpectralTerms(
-      refFreq, isLogarithmic,
-      xt::concatenate(xt::xtuple(spectrum_.GetSpectralTerms(), new_terms), 0));
-}
 
+  spectrum_.SetSpectralTerms(refFreq, isLogarithmic,
+                             xt::adapt(std::vector<double>(first, last)));
+};
 } // namespace base
 } // namespace dp3
 
diff --git a/include/PointSourceCollection.h b/include/PointSourceCollection.h
index e4cfebf2ac8148b2c208677b5cf97e01ec190c1d..8171a83332d2daabc901dd34d86b42d6304c16d8 100644
--- a/include/PointSourceCollection.h
+++ b/include/PointSourceCollection.h
@@ -22,15 +22,12 @@ using PointSource = dp3::base::PointSource;
 class PointSourceCollection : public ObjectCollection<PointSource> {
 public:
   Directions direction_vector;
-  std::vector<StokesVector> spectrum_vector;
   std::vector<Spectrum> spectrums;
   std::vector<size_t> beam_id;
   std::set<size_t> beam_ids;
-  StokesVector stokes_vector;
 
   void Add(const PointSource &point_source) {
     direction_vector.Add(point_source.GetDirection());
-    stokes_vector.Add(point_source.GetStokes());
     spectrums.push_back(point_source.GetSpectrum());
     beam_id.push_back(point_source.GetBeamId());
     beam_ids.insert(point_source.GetBeamId());
@@ -38,20 +35,32 @@ public:
 
   void Reserve(size_t new_size) {
     direction_vector.Reserve(new_size);
-    stokes_vector.Reserve(new_size);
     spectrums.reserve(new_size);
     beam_id.reserve(new_size);
   }
 
   void Clear() {
-    PointSourceCollection::Clear();
     direction_vector.Clear();
-    spectrum_vector.clear();
     spectrums.clear();
     beam_id.clear();
     beam_ids.clear();
   }
 
+  PointSource operator[](size_t i) const {
+    return PointSource(direction_vector[i], spectrums[i], beam_id[i]);
+  }
+
+  std::unique_ptr<PointSourceCollection> SelectBeamID(size_t beam_id) {
+    auto selected = std::make_unique<PointSourceCollection>();
+
+    for (size_t i = 0; i < Size(); ++i) {
+      if (this->beam_id[i] == beam_id) {
+        selected->Add(operator[](i));
+      }
+    }
+    return std::move(selected);
+  }
+
   size_t Size() const { return direction_vector.Size(); }
 };
 
diff --git a/include/Predict.h b/include/Predict.h
index 0851580c2d1f822e53311ea108a0f52098e5f0ce..08f2579d543d5b75ad46b35b115c9b9546cf53fa 100644
--- a/include/Predict.h
+++ b/include/Predict.h
@@ -86,7 +86,7 @@ inline void spectrum(const PointSource &component, size_t nChannel,
 /// @{
 
 /**
- * @brief Simulator to compute visibilities given a sky model
+ * @brief Predict class to compute visibilities given a sky model
  *
  * This class computes visibilities given model components in the sky.
  * Effectively, it evaluates the equation:
@@ -126,6 +126,11 @@ public:
            xt::xtensor<dcomplex, 3> &buffer) const;
   void run(PredictPlanExec &plan, const GaussianSourceCollection &sources,
            xt::xtensor<dcomplex, 3> &buffer) const;
+
+  void run(PredictPlanExec &plan, const GaussianSourceCollection &sources,
+           xt::xtensor<dcomplex, 4> &buffer) const;
+  void run(PredictPlanExec &plan, const PointSourceCollection &sources,
+           xt::xtensor<dcomplex, 4> &buffer) const;
 };
 
 /// @}
diff --git a/include/Spectrum.h b/include/Spectrum.h
index 3fabc95f515cd20dc63e1916640ae4e0c168fdf5..369a4ddf14a8e62e01f6031ea08a9a87a090f697 100644
--- a/include/Spectrum.h
+++ b/include/Spectrum.h
@@ -11,7 +11,6 @@
 
 #include <Stokes.h>
 #include <StokesVector.h>
-
 #include <casacore/casa/BasicSL/Constants.h>
 #include <complex>
 #include <xtensor/xcomplex.hpp>
@@ -91,7 +90,7 @@ public:
     return has_logarithmic_spectral_index_;
   }
   bool HasSpectralTerms() const { return spectral_terms_.size() > 0; }
-
+  bool HasRotationMeasure() const { return has_rotation_measure_; }
   void SetPolarization(double angle, double factor) {
     polarization_angle_ = angle;
     polarization_factor_ = factor;
@@ -207,7 +206,6 @@ private:
     return EvaluatePolynomial(x, parameters) * x;
   }
 
-private:
   Stokes reference_flux_;
 
   xt::xtensor<double, 1> spectral_terms_;
diff --git a/src/GaussianSource.cpp b/src/GaussianSource.cpp
index 2e558aaaf947679d0161742d26d0d80624943964..3d67e4034913e3ea60dbd8db2ccd0c72f71795ab 100644
--- a/src/GaussianSource.cpp
+++ b/src/GaussianSource.cpp
@@ -20,6 +20,16 @@ GaussianSource::GaussianSource(const Direction &direction, const Stokes &stokes,
     : PointSource(direction, stokes, beam_id), position_angle_(0.0),
       is_position_angle_absolute_(true), major_axis_(0.0), minor_axis_(0.0) {}
 
+GaussianSource::GaussianSource(const Direction &direction,
+                               const Spectrum &spectrum, double position_angle,
+                               bool is_position_angle_absolute,
+                               double minor_axis, double major_axis,
+                               size_t beam_id)
+    : PointSource(direction, spectrum, beam_id),
+      position_angle_(position_angle),
+      is_position_angle_absolute_(is_position_angle_absolute),
+      major_axis_(major_axis), minor_axis_(minor_axis) {}
+
 void GaussianSource::SetPositionAngle(double angle) { position_angle_ = angle; }
 
 void GaussianSource::SetMajorAxis(double fwhm) { major_axis_ = fwhm; }
diff --git a/src/PointSource.cpp b/src/PointSource.cpp
index 569edda8ae3e1fd59af3e3275a8a113d9fc2b99e..88c1ef96dc60f6fc2e2feef6426694c447d5638b 100644
--- a/src/PointSource.cpp
+++ b/src/PointSource.cpp
@@ -18,13 +18,14 @@ namespace base {
 
 PointSource::PointSource(const Direction &position, const Stokes &stokes,
                          const size_t beam_id)
-    : direction_(position), stokes_(stokes), reference_frequency_(0.0),
-      polarizated_fraction_(0.0), polarization_angle_(0.0),
-      rotation_measure_(0.0), has_rotation_measure_(false),
-      has_logarithmic_si_(true), beam_id_(beam_id) {
+    : direction_(position), beam_id_(beam_id) {
   spectrum_.SetReferenceFlux(stokes);
 }
 
+PointSource::PointSource(const Direction &position, const Spectrum &spectrum,
+                         const size_t beam_id)
+    : direction_(position), beam_id_(beam_id), spectrum_(spectrum) {}
+
 void PointSource::SetDirection(const Direction &direction) {
   direction_ = direction;
 }
@@ -36,31 +37,30 @@ void PointSource::ComputeSpectrum(
 }
 
 void PointSource::SetRotationMeasure(double fraction, double angle, double rm) {
-  polarizated_fraction_ = fraction;
-  polarization_angle_ = angle;
-  rotation_measure_ = rm;
-  has_rotation_measure_ = true;
+  spectrum_.SetPolarization(angle, fraction);
+  spectrum_.SetRotationMeasure(rm);
 }
 
 // FIXME: legacy code, remove it later
 Stokes PointSource::GetStokes(double freq) const {
-  Stokes stokes(stokes_);
+  Stokes stokes(spectrum_.GetReferenceFlux());
 
   if (HasSpectralTerms()) {
-    if (has_logarithmic_si_) {
+    if (spectrum_.HasLogarithmicSpectralIndex()) {
       // Compute spectral index as:
       // (v / v0) ^ (c0 + c1 * log10(v / v0) + c2 * log10(v / v0)^2 + ...)
       // Where v is the frequency and v0 is the reference frequency.
 
       // Compute log10(v / v0).
-      double base = log10(freq) - log10(reference_frequency_);
+      double base = log10(freq) - log10(spectrum_.GetReferenceFrequency());
 
       // Compute c0 + log10(v / v0) * c1 + log10(v / v0)^2 * c2 + ...
       // using Horner's rule.
       double exponent = 0.0;
-      typedef std::vector<double>::const_reverse_iterator iterator_type;
-      for (iterator_type it = spectral_terms_.rbegin(),
-                         end = spectral_terms_.rend();
+      typedef xt::xtensor<double, 1>::const_reverse_iterator iterator_type;
+      auto spectral_terms = spectrum_.GetSpectralTerms();
+      for (iterator_type it = spectral_terms.rbegin(),
+                         end = spectral_terms.rend();
            it != end; ++it) {
         exponent = exponent * base + *it;
       }
@@ -70,9 +70,10 @@ Stokes PointSource::GetStokes(double freq) const {
       stokes.I *= pow(10., base * exponent);
       stokes.V *= pow(10., base * exponent);
     } else {
-      double x = freq / reference_frequency_ - 1.0;
-      typedef std::vector<double>::const_reverse_iterator iterator_type;
+      double x = freq / spectrum_.GetReferenceFrequency() - 1.0;
+      typedef xt::xtensor<double, 1>::const_reverse_iterator iterator_type;
       double val = 0.0;
+      auto spectral_terms_ = spectrum_.GetSpectralTerms();
       for (iterator_type it = spectral_terms_.rbegin(),
                          end = spectral_terms_.rend();
            it != end; ++it) {
@@ -85,18 +86,22 @@ Stokes PointSource::GetStokes(double freq) const {
 
   if (HasRotationMeasure()) {
     double lambda = casacore::C::c / freq;
-    double chi =
-        2.0 * (polarization_angle_ + rotation_measure_ * lambda * lambda);
-    double stokesQU = stokes.I * polarizated_fraction_;
+    double chi = 2.0 * (spectrum_.GetPolarizationAngle() +
+                        spectrum_.GetRotationMeasure() * lambda * lambda);
+    double stokesQU = stokes.I * spectrum_.GetPolarizationFactor();
     stokes.Q = stokesQU * cos(chi);
     stokes.U = stokesQU * sin(chi);
   }
   return stokes;
 }
 
-bool PointSource::HasSpectralTerms() const { return !spectral_terms_.empty(); }
+bool PointSource::HasSpectralTerms() const {
+  return spectrum_.HasSpectralTerms();
+}
 
-bool PointSource::HasRotationMeasure() const { return has_rotation_measure_; }
+bool PointSource::HasRotationMeasure() const {
+  return spectrum_.HasRotationMeasure();
+}
 
 } // namespace base
 } // namespace dp3
diff --git a/src/Predict.cpp b/src/Predict.cpp
index a05045a49f4583a27e98bba51fb83a1875b7835a..1ed984de0413e9d57f49e731b6cef537d05935f0 100644
--- a/src/Predict.cpp
+++ b/src/Predict.cpp
@@ -34,5 +34,40 @@ void Predict::run(PredictPlanExec &plan,
   plan.Compute(sources, buffer);
 }
 
+void Predict::run(PredictPlanExec &plan,
+                  const GaussianSourceCollection &sources,
+                  xt::xtensor<dcomplex, 4> &buffer) const {
+  plan.Precompute(sources);
+
+  xt::xtensor<dcomplex, 3>::shape_type tmp_buffer_shape = {
+      buffer.shape(1), buffer.shape(2), buffer.shape(3)};
+  for (auto &beam_id : sources.beam_ids) {
+    xt::xtensor<dcomplex, 3> direction_buffer(tmp_buffer_shape);
+
+    // TODO compute should be extended to accumulate per beam_id
+    plan.Compute(sources, direction_buffer);
+    xt::view(buffer, beam_id, xt::all(), xt::all(), xt::all()) =
+        direction_buffer;
+  }
+}
+
+void Predict::run(PredictPlanExec &plan, const PointSourceCollection &sources,
+                  xt::xtensor<dcomplex, 4> &buffer) const {
+  // buffer dimensions are (nBean, nFreq, nBaselines, nSources)
+  plan.Precompute(sources);
+
+  xt::xtensor<dcomplex, 3>::shape_type tmp_buffer_shape = {
+      buffer.shape(1), buffer.shape(2), buffer.shape(3)};
+  for (auto &beam_id : sources.beam_ids) {
+    xt::xtensor<dcomplex, 3> direction_buffer(tmp_buffer_shape);
+
+    // TODO compute should be extended to accumulate per beam_id
+    plan.Compute(sources, direction_buffer);
+    xt::view(buffer, beam_id, xt::all(), xt::all(), xt::all()) =
+        direction_buffer;
+  }
+}
+// plan.Compute(sources, buffer);}
+
 } // namespace base
 } // namespace dp3
diff --git a/tests/test_gaussiansourcecollection.cpp b/tests/test_gaussiansourcecollection.cpp
index e797036c0ff66bbf2f11eaf254454401ccd77f95..3749a1c31c1b1ec487d68b9ff59a5b71db08aaa2 100644
--- a/tests/test_gaussiansourcecollection.cpp
+++ b/tests/test_gaussiansourcecollection.cpp
@@ -14,10 +14,10 @@ BOOST_AUTO_TEST_CASE(add_single) {
   collection.Add(source);
   BOOST_CHECK_EQUAL(collection.Size(), 1);
   BOOST_CHECK_EQUAL(collection.direction_vector.Size(), 1);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.Size(), 1);
+
   BOOST_CHECK_EQUAL(collection.direction_vector.ra[0], 0.1);
   BOOST_CHECK_EQUAL(collection.direction_vector.dec[0], 0.2);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.I[0], 1.0);
+  BOOST_CHECK_EQUAL(collection.spectrums[0].GetReferenceFlux().I, 1.0);
 }
 
 BOOST_AUTO_TEST_CASE(reserve) {
@@ -27,13 +27,13 @@ BOOST_AUTO_TEST_CASE(reserve) {
   PointSourceCollection collection;
   collection.Reserve(2);
 
-  auto ptr = collection.stokes_vector.I.data();
+  auto ptr = collection.direction_vector.ra.data();
 
   collection.Add(sources[0]);
   collection.Add(sources[1]);
 
   BOOST_CHECK_EQUAL(collection.Size(), 2);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.I.data(), ptr);
+  BOOST_CHECK_EQUAL(collection.direction_vector.ra.data(), ptr);
 }
 
 BOOST_AUTO_TEST_SUITE_END()
\ No newline at end of file
diff --git a/tests/test_pointsourcecollection.cpp b/tests/test_pointsourcecollection.cpp
index dcddb1e7014d9bc73e4d4f8e9b67c5e439a0f7b9..2d48acc4689a8468c29030da4945725e3934f3af 100644
--- a/tests/test_pointsourcecollection.cpp
+++ b/tests/test_pointsourcecollection.cpp
@@ -14,10 +14,10 @@ BOOST_AUTO_TEST_CASE(add_single) {
   collection.Add(source);
   BOOST_CHECK_EQUAL(collection.Size(), 1);
   BOOST_CHECK_EQUAL(collection.direction_vector.Size(), 1);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.Size(), 1);
+  BOOST_CHECK_EQUAL(collection.spectrums.size(), 1);
   BOOST_CHECK_EQUAL(collection.direction_vector.ra[0], 0.1);
   BOOST_CHECK_EQUAL(collection.direction_vector.dec[0], 0.2);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.I[0], 1.0);
+  BOOST_CHECK_EQUAL(collection.spectrums[0].GetReferenceFlux().I, 1.0);
 }
 
 BOOST_AUTO_TEST_CASE(add_single_unspecified_beam_id) {
@@ -26,10 +26,10 @@ BOOST_AUTO_TEST_CASE(add_single_unspecified_beam_id) {
   collection.Add(source);
   BOOST_CHECK_EQUAL(collection.Size(), 1);
   BOOST_CHECK_EQUAL(collection.direction_vector.Size(), 1);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.Size(), 1);
+  BOOST_CHECK_EQUAL(collection.spectrums.size(), 1);
   BOOST_CHECK_EQUAL(collection.direction_vector.ra[0], 0.1);
   BOOST_CHECK_EQUAL(collection.direction_vector.dec[0], 0.2);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.I[0], 1.0);
+
   BOOST_CHECK_EQUAL(collection.beam_id[0], 0);
 }
 
@@ -39,10 +39,10 @@ BOOST_AUTO_TEST_CASE(add_single_with_beam_id) {
   collection.Add(source);
   BOOST_CHECK_EQUAL(collection.Size(), 1);
   BOOST_CHECK_EQUAL(collection.direction_vector.Size(), 1);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.Size(), 1);
+  BOOST_CHECK_EQUAL(collection.spectrums.size(), 1);
   BOOST_CHECK_EQUAL(collection.direction_vector.ra[0], 0.1);
   BOOST_CHECK_EQUAL(collection.direction_vector.dec[0], 0.2);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.I[0], 1.0);
+  BOOST_CHECK_EQUAL(collection.spectrums[0].GetReferenceFlux().I, 1.0);
   BOOST_CHECK_EQUAL(collection.beam_id[0], 15);
   BOOST_CHECK_EQUAL(collection.beam_ids.size(), 1);
 }
@@ -69,13 +69,13 @@ BOOST_AUTO_TEST_CASE(reserve) {
   PointSourceCollection collection;
   collection.Reserve(2);
 
-  auto ptr = collection.stokes_vector.I.data();
+  auto ptr = collection.direction_vector.ra.data();
 
   collection.Add(sources[0]);
   collection.Add(sources[1]);
 
   BOOST_CHECK_EQUAL(collection.Size(), 2);
-  BOOST_CHECK_EQUAL(collection.stokes_vector.I.data(), ptr);
+  BOOST_CHECK_EQUAL(collection.direction_vector.ra.data(), ptr);
 }
 
 BOOST_AUTO_TEST_SUITE_END()
\ No newline at end of file
diff --git a/tests/test_spectrum.cpp b/tests/test_spectrum.cpp
index 10badfb3df0f39e90b90baf4ee09876255802d3a..2131121a89a3b03f8c6cb76490d69531808acd27 100644
--- a/tests/test_spectrum.cpp
+++ b/tests/test_spectrum.cpp
@@ -24,13 +24,11 @@ struct SpectrumFixture {
 };
 
 BOOST_FIXTURE_TEST_CASE(normal_no_rotation, SpectrumFixture) {
-  dp3::base::PointSource point_source(dp3::base::Direction(0.0, 0.0),
-                                      spectrum.GetReferenceFlux());
+  dp3::base::PointSource point_source(dp3::base::Direction(0.0, 0.0), spectrum);
 
   point_source.SetSpectralTerms(
       spectrum.GetReferenceFrequency(), spectrum.HasLogarithmicSpectralIndex(),
       spectrum.GetSpectralTerms().begin(), spectrum.GetSpectralTerms().end());
-  point_source.SetRotationMeasure(0.0, 0.0, 0.0);
 
   Stokes stokes_nu0(0, 0, 0, 0);
   Stokes stokes_nu1(0, 0, 0, 0);