diff --git a/conanfile.txt b/conanfile.txt
index f7e1f43711091dfeaec269ed86395720bd4553cc..348ff6e5309bd20e68dd2edfefe3454228da5a58 100644
--- a/conanfile.txt
+++ b/conanfile.txt
@@ -4,7 +4,7 @@ catch2/2.13.8
 eigen/3.3.8
 boost/1.78.0
 zlib/1.2.13
-proposal/7.6.1
+proposal/7.6.2
 yaml-cpp/0.7.0
 arrow/10.0.0
 cli11/1.9.1
diff --git a/corsika/detail/modules/proposal/ContinuousProcess.inl b/corsika/detail/modules/proposal/ContinuousProcess.inl
index d0accf113dea87805fb0f3598046faab838e3550..b190619b39b88f0b81ceafbb7fff4a73992d4862 100644
--- a/corsika/detail/modules/proposal/ContinuousProcess.inl
+++ b/corsika/detail/modules/proposal/ContinuousProcess.inl
@@ -34,18 +34,15 @@ namespace corsika::proposal {
         code); //! energy resolutions globally defined for individual particles
     auto c = p_cross->second(media.at(comp.getHash()), emCut);
 
-    // Use higland multiple scattering and deactivate stochastic deflection by
-    // passing an empty vector
+    // choose multiple scattering model
     static constexpr auto ms_type = PROPOSAL::MultipleScatteringType::MoliereInterpol;
-    auto s_type = std::vector<PROPOSAL::InteractionType>();
 
     // Build displacement integral and scattering object and interpolate them too and
     // saved in the calc map by a key build out of a hash of composed of the component and
     // particle code.
-    auto calculator =
-        Calculator{PROPOSAL::make_displacement(c, true),
-                   PROPOSAL::make_scattering(ms_type, s_type, particle[code],
-                                             media.at(comp.getHash()))};
+    auto calculator = Calculator{PROPOSAL::make_displacement(c, true),
+                                 PROPOSAL::make_multiple_scattering(
+                                     ms_type, particle[code], media.at(comp.getHash()))};
     calc[std::make_pair(comp.getHash(), code)] = std::move(calculator);
   }
 
@@ -65,33 +62,47 @@ namespace corsika::proposal {
     // get or build corresponding calculators
     auto c = getCalculator(step.getParticlePre(), calc);
 
-    // Cast corsika vector to proposal vector
-    auto particle_dir = step.getDirectionPre();
-    auto d = particle_dir.getComponents();
-    auto direction = PROPOSAL::Cartesian3D(d.getX().magnitude(), d.getY().magnitude(),
-                                           d.getZ().magnitude());
+    auto initial_particle_dir = step.getDirectionPre();
 
     auto E_i_total = step.getEkinPre() + step.getParticlePre().getMass();
     auto E_f_total = E_i_total - loss;
 
-    // draw random numbers required for scattering process
+    // sample scattering_angle, which is the combination sqrt(theta_x^2 + theta_y^2),
+    // where theta_x and theta_y itself are independent angles drawn from the multiple
+    // scattering distribution
     std::uniform_real_distribution<double> distr(0., 1.);
-    auto rnd = std::array<double, 4>();
-    for (auto& it : rnd) it = distr(RNG_);
+    auto scattering_angle = (c->second).scatter->CalculateScatteringAngle2D(
+        grammage / 1_g * square(1_cm), E_i_total / 1_MeV, E_f_total / 1_MeV, distr(RNG_),
+        distr(RNG_));
+
+    auto const& root = initial_particle_dir.getCoordinateSystem();
+
+    // construct vector that is normal to initial direction.
+    DirectionVector normal_vec{root, {0, 0, 0}};
+    if (initial_particle_dir.getX(root) > 0.1) {
+      normal_vec = {
+          root, {-initial_particle_dir.getY(root), initial_particle_dir.getX(root), 0}};
+    } else {
+      // if x is small, use y and z to construct normal vector
+      normal_vec = {
+          root, {0, -initial_particle_dir.getZ(root), initial_particle_dir.getY(root)}};
+    }
+
+    // rotation of zenith by moliere_angle
+    CoordinateSystemPtr rotated1 =
+        make_rotation(root, normal_vec.getComponents(), scattering_angle);
 
-    // calculate deflection based on particle energy, loss
-    auto deflection = (c->second).scatter->CalculateMultipleScattering(
-        grammage / 1_g * square(1_cm), E_i_total / 1_MeV, E_f_total / 1_MeV, rnd);
+    // rotation of azimuth by random angle between 0 and 2*PI
+    std::uniform_real_distribution<double> distr_azimuth(0., 2 * M_PI);
+    CoordinateSystemPtr rotated2 = make_rotation(
+        rotated1, initial_particle_dir.getComponents(), distr_azimuth(RNG_));
 
-    [[maybe_unused]] auto [unused1, final_direction] =
-        PROPOSAL::multiple_scattering::ScatterInitialDirection(direction, deflection);
+    DirectionVector scattered_particle_dir{root,
+                                           initial_particle_dir.getComponents(rotated2)};
 
     // update particle direction after continuous loss caused by multiple
     // scattering
-    DirectionVector dU_{
-        particle_dir.getCoordinateSystem(),
-        {final_direction.GetX(), final_direction.GetY(), final_direction.GetZ()}};
-    DirectionVector diff_dir_ = dU_ - particle_dir;
+    DirectionVector diff_dir_ = scattered_particle_dir - initial_particle_dir;
     step.add_dU(diff_dir_);
   }
 
diff --git a/corsika/modules/proposal/ContinuousProcess.hpp b/corsika/modules/proposal/ContinuousProcess.hpp
index 054df34d34e2f4e7fa6526f04ec89da504f2fbf9..fc597be321f100bf6348babb718d03df326f7ef0 100644
--- a/corsika/modules/proposal/ContinuousProcess.hpp
+++ b/corsika/modules/proposal/ContinuousProcess.hpp
@@ -39,7 +39,7 @@ namespace corsika::proposal {
 
     struct Calculator {
       std::unique_ptr<PROPOSAL::Displacement> disp;
-      std::unique_ptr<PROPOSAL::Scattering> scatter;
+      std::unique_ptr<PROPOSAL::multiple_scattering::Parametrization> scatter;
     };
 
     std::unordered_map<calc_key_t, Calculator, hash>