diff --git a/crates/bevy_animation/src/lib.rs b/crates/bevy_animation/src/lib.rs
index 02f8139a697d5..a389da5176d32 100644
--- a/crates/bevy_animation/src/lib.rs
+++ b/crates/bevy_animation/src/lib.rs
@@ -2,7 +2,7 @@
 
 #![warn(missing_docs)]
 
-use std::ops::Deref;
+use std::ops::{Add, Deref, Mul};
 use std::time::Duration;
 
 use bevy_app::{App, Plugin, PostUpdate};
@@ -126,7 +126,7 @@ impl VariableCurve {
 }
 
 /// Interpolation method to use between keyframes.
-#[derive(Reflect, Clone, Debug)]
+#[derive(Reflect, Copy, Clone, Debug)]
 pub enum Interpolation {
     /// Linear interpolation between the two closest keyframes.
     Linear,
@@ -137,6 +137,69 @@ pub enum Interpolation {
     CubicSpline,
 }
 
+impl Interpolation {
+    #[inline(always)]
+    fn interpolate<T: Interpolatable>(
+        self,
+        lerp: f32,
+        duration: f32,
+        keyframes: &[T],
+        stride: usize,
+    ) -> T {
+        match self {
+            Interpolation::Linear => {
+                let value_start = keyframes[stride * 0];
+                let value_end = keyframes[stride * 1];
+                value_start._lerp(value_end, lerp)
+            }
+            Interpolation::Step => keyframes[stride * 0],
+            Interpolation::CubicSpline => {
+                let value_start = keyframes[stride * 1];
+                let tangent_out_start = keyframes[stride * 2];
+                let tangent_in_end = keyframes[stride * 3];
+                let value_end = keyframes[stride * 4];
+                (value_start * (2.0 * lerp.powi(3) - 3.0 * lerp.powi(2) + 1.0))
+                    + (tangent_out_start * (duration * (lerp.powi(3) - 2.0 * lerp.powi(2) + lerp)))
+                    + (value_end * (-2.0 * lerp.powi(3) + 3.0 * lerp.powi(2)))
+                    + (tangent_in_end * (duration * (lerp.powi(3) - lerp.powi(2))))
+            }
+        }
+    }
+
+    #[inline(always)]
+    fn step_size(self) -> usize {
+        match self {
+            Interpolation::Linear | Interpolation::Step => 1,
+            Interpolation::CubicSpline => 3,
+        }
+    }
+}
+
+trait Interpolatable: Sized + Copy + Mul<f32, Output = Self> + Add<Self, Output = Self> {
+    fn _lerp(self, rhs: Self, t: f32) -> Self;
+}
+
+impl Interpolatable for f32 {
+    #[inline(always)]
+    fn _lerp(self, rhs: Self, t: f32) -> Self {
+        self.lerp(rhs, t)
+    }
+}
+
+impl Interpolatable for Vec3 {
+    #[inline(always)]
+    fn _lerp(self, rhs: Self, t: f32) -> Self {
+        self.lerp(rhs, t)
+    }
+}
+
+impl Interpolatable for Quat {
+    #[inline(always)]
+    fn _lerp(self, rhs: Self, t: f32) -> Self {
+        self.slerp(rhs, t)
+    }
+}
+
 /// Path to an entity, with [`Name`]s. Each entity in a path must have a name.
 #[derive(Reflect, Clone, Debug, Hash, PartialEq, Eq, Default)]
 pub struct EntityPath {
@@ -648,42 +711,6 @@ fn run_animation_player(
     }
 }
 
-/// Update `weights` based on weights in `keyframe` with a linear interpolation
-/// on `key_lerp`.
-fn lerp_morph_weights(weights: &mut [f32], keyframe: impl Iterator<Item = f32>, key_lerp: f32) {
-    let zipped = weights.iter_mut().zip(keyframe);
-    for (morph_weight, keyframe) in zipped {
-        *morph_weight = morph_weight.lerp(keyframe, key_lerp);
-    }
-}
-
-/// Extract a keyframe from a list of keyframes by index.
-///
-/// # Panics
-///
-/// When `key_index * target_count` is larger than `keyframes`
-///
-/// This happens when `keyframes` is not formatted as described in
-/// [`Keyframes::Weights`]. A possible cause is [`AnimationClip`] not being
-/// meant to be used for the [`MorphWeights`] of the entity it's being applied to.
-fn get_keyframe(target_count: usize, keyframes: &[f32], key_index: usize) -> &[f32] {
-    let start = target_count * key_index;
-    let end = target_count * (key_index + 1);
-    &keyframes[start..end]
-}
-
-// Helper macro for cubic spline interpolation
-// it needs to work on `f32`, `Vec3` and `Quat`
-// TODO: replace by a function if the proper trait bounds can be figured out
-macro_rules! cubic_spline_interpolation {
-    ($value_start: expr, $tangent_out_start: expr, $tangent_in_end: expr, $value_end: expr, $lerp: expr, $step_duration: expr,) => {
-        $value_start * (2.0 * $lerp.powi(3) - 3.0 * $lerp.powi(2) + 1.0)
-            + $tangent_out_start * ($step_duration) * ($lerp.powi(3) - 2.0 * $lerp.powi(2) + $lerp)
-            + $value_end * (-2.0 * $lerp.powi(3) + 3.0 * $lerp.powi(2))
-            + $tangent_in_end * ($step_duration) * ($lerp.powi(3) - $lerp.powi(2))
-    };
-}
-
 #[allow(clippy::too_many_arguments)]
 fn apply_animation(
     weight: f32,
@@ -760,11 +787,13 @@ fn apply_animation(
                     Keyframes::Weights(keyframes) => {
                         if let Some(morphs) = &mut morphs {
                             let target_count = morphs.weights().len();
-                            lerp_morph_weights(
-                                morphs.weights_mut(),
-                                get_keyframe(target_count, keyframes, 0).iter().copied(),
-                                weight,
-                            );
+                            for (x, y) in morphs
+                                .weights_mut()
+                                .iter_mut()
+                                .zip(&keyframes[..target_count])
+                            {
+                                *x = x.lerp(*y, weight)
+                            }
                         }
                     }
                 }
@@ -808,128 +837,57 @@ fn apply_keyframe(
     transform: &mut Mut<Transform>,
     morphs: &mut Option<Mut<MorphWeights>>,
 ) {
-    match (&curve.interpolation, &curve.keyframes) {
-        (Interpolation::Step, Keyframes::Rotation(keyframes)) => {
-            transform.rotation = transform.rotation.slerp(keyframes[step_start], weight);
-        }
-        (Interpolation::Linear, Keyframes::Rotation(keyframes)) => {
-            let rot_start = keyframes[step_start];
-            let mut rot_end = keyframes[step_start + 1];
-            // Choose the smallest angle for the rotation
-            if rot_end.dot(rot_start) < 0.0 {
-                rot_end = -rot_end;
-            }
-            // Rotations are using a spherical linear interpolation
-            let rot = rot_start.normalize().slerp(rot_end.normalize(), lerp);
-            transform.rotation = transform.rotation.slerp(rot, weight);
-        }
-        (Interpolation::CubicSpline, Keyframes::Rotation(keyframes)) => {
-            let value_start = keyframes[step_start * 3 + 1];
-            let tangent_out_start = keyframes[step_start * 3 + 2];
-            let tangent_in_end = keyframes[(step_start + 1) * 3];
-            let value_end = keyframes[(step_start + 1) * 3 + 1];
-            let result = cubic_spline_interpolation!(
-                value_start,
-                tangent_out_start,
-                tangent_in_end,
-                value_end,
-                lerp,
-                duration,
+    match &curve.keyframes {
+        Keyframes::Rotation(keyframes) => {
+            transform.rotation = transform.rotation.slerp(
+                curve.interpolation.interpolate(
+                    lerp,
+                    duration,
+                    &keyframes[step_start * curve.interpolation.step_size()..],
+                    1,
+                ),
+                weight,
             );
-            transform.rotation = transform.rotation.slerp(result.normalize(), weight);
         }
-        (Interpolation::Step, Keyframes::Translation(keyframes)) => {
-            transform.translation = transform.translation.lerp(keyframes[step_start], weight);
-        }
-        (Interpolation::Linear, Keyframes::Translation(keyframes)) => {
-            let translation_start = keyframes[step_start];
-            let translation_end = keyframes[step_start + 1];
-            let result = translation_start.lerp(translation_end, lerp);
-            transform.translation = transform.translation.lerp(result, weight);
-        }
-        (Interpolation::CubicSpline, Keyframes::Translation(keyframes)) => {
-            let value_start = keyframes[step_start * 3 + 1];
-            let tangent_out_start = keyframes[step_start * 3 + 2];
-            let tangent_in_end = keyframes[(step_start + 1) * 3];
-            let value_end = keyframes[(step_start + 1) * 3 + 1];
-            let result = cubic_spline_interpolation!(
-                value_start,
-                tangent_out_start,
-                tangent_in_end,
-                value_end,
-                lerp,
-                duration,
+        Keyframes::Translation(keyframes) => {
+            transform.translation = transform.translation.lerp(
+                curve.interpolation.interpolate(
+                    lerp,
+                    duration,
+                    &keyframes[step_start * curve.interpolation.step_size()..],
+                    1,
+                ),
+                weight,
             );
-            transform.translation = transform.translation.lerp(result, weight);
-        }
-        (Interpolation::Step, Keyframes::Scale(keyframes)) => {
-            transform.scale = transform.scale.lerp(keyframes[step_start], weight);
         }
-        (Interpolation::Linear, Keyframes::Scale(keyframes)) => {
-            let scale_start = keyframes[step_start];
-            let scale_end = keyframes[step_start + 1];
-            let result = scale_start.lerp(scale_end, lerp);
-            transform.scale = transform.scale.lerp(result, weight);
-        }
-        (Interpolation::CubicSpline, Keyframes::Scale(keyframes)) => {
-            let value_start = keyframes[step_start * 3 + 1];
-            let tangent_out_start = keyframes[step_start * 3 + 2];
-            let tangent_in_end = keyframes[(step_start + 1) * 3];
-            let value_end = keyframes[(step_start + 1) * 3 + 1];
-            let result = cubic_spline_interpolation!(
-                value_start,
-                tangent_out_start,
-                tangent_in_end,
-                value_end,
-                lerp,
-                duration,
+        Keyframes::Scale(keyframes) => {
+            transform.scale = transform.scale.lerp(
+                curve.interpolation.interpolate(
+                    lerp,
+                    duration,
+                    &keyframes[step_start * curve.interpolation.step_size()..],
+                    1,
+                ),
+                weight,
             );
-            transform.scale = transform.scale.lerp(result, weight);
-        }
-        (Interpolation::Step, Keyframes::Weights(keyframes)) => {
-            if let Some(morphs) = morphs {
-                let target_count = morphs.weights().len();
-                let morph_start = get_keyframe(target_count, keyframes, step_start);
-                lerp_morph_weights(morphs.weights_mut(), morph_start.iter().copied(), weight);
-            }
         }
-        (Interpolation::Linear, Keyframes::Weights(keyframes)) => {
+        Keyframes::Weights(keyframes) => {
             if let Some(morphs) = morphs {
                 let target_count = morphs.weights().len();
-                let morph_start = get_keyframe(target_count, keyframes, step_start);
-                let morph_end = get_keyframe(target_count, keyframes, step_start + 1);
-                let result = morph_start
-                    .iter()
-                    .zip(morph_end)
-                    .map(|(a, b)| a.lerp(*b, lerp));
-                lerp_morph_weights(morphs.weights_mut(), result, weight);
-            }
-        }
-        (Interpolation::CubicSpline, Keyframes::Weights(keyframes)) => {
-            if let Some(morphs) = morphs {
-                let target_count = morphs.weights().len();
-                let morph_start = get_keyframe(target_count, keyframes, step_start * 3 + 1);
-                let tangents_out_start = get_keyframe(target_count, keyframes, step_start * 3 + 2);
-                let tangents_in_end = get_keyframe(target_count, keyframes, (step_start + 1) * 3);
-                let morph_end = get_keyframe(target_count, keyframes, (step_start + 1) * 3 + 1);
-                let result = morph_start
-                    .iter()
-                    .zip(tangents_out_start)
-                    .zip(tangents_in_end)
-                    .zip(morph_end)
-                    .map(
-                        |(((value_start, tangent_out_start), tangent_in_end), value_end)| {
-                            cubic_spline_interpolation!(
-                                value_start,
-                                tangent_out_start,
-                                tangent_in_end,
-                                value_end,
-                                lerp,
-                                duration,
-                            )
-                        },
+                for (morph, morph_weight) in morphs.weights_mut().iter_mut().enumerate() {
+                    *morph_weight = morph_weight.lerp(
+                        curve.interpolation.interpolate(
+                            lerp,
+                            duration,
+                            &keyframes[step_start
+                                * curve.interpolation.step_size()
+                                * target_count
+                                + morph..],
+                            target_count,
+                        ),
+                        weight,
                     );
-                lerp_morph_weights(morphs.weights_mut(), result, weight);
+                }
             }
         }
     }