Animation: Butterworth Smoothing filter #106952

Merged
Christoph Lendenfeld merged 56 commits from ChrisLend/blender:graph_editor_butterworth into main 2023-07-13 09:10:49 +02:00
6 changed files with 512 additions and 7 deletions

View File

@ -309,6 +309,7 @@ class GRAPH_MT_key_smoothing(Menu):
layout.operator_context = 'INVOKE_DEFAULT'
layout.operator("graph.gaussian_smooth", text="Smooth (Gaussian)")
layout.operator("graph.smooth", text="Smooth (Legacy)")
layout.operator("graph.butterworth_smooth")
class GRAPH_MT_key(Menu):

View File

@ -387,6 +387,197 @@ void blend_to_default_fcurve(PointerRNA *id_ptr, FCurve *fcu, const float factor
BKE_fcurve_keyframe_move_value_with_handles(&fcu->bezt[i], key_y_value);
}
}
/* ---------------- */
ButterworthCoefficients *ED_anim_allocate_butterworth_coefficients(const int filter_order)
{
ButterworthCoefficients *bw_coeff = static_cast<ButterworthCoefficients *>(
MEM_callocN(sizeof(ButterworthCoefficients), "Butterworth Coefficients"));
bw_coeff->filter_order = filter_order;
bw_coeff->d1 = static_cast<double *>(
MEM_callocN(sizeof(double) * filter_order, "coeff filtered"));
bw_coeff->d2 = static_cast<double *>(
MEM_callocN(sizeof(double) * filter_order, "coeff samples"));
bw_coeff->A = static_cast<double *>(MEM_callocN(sizeof(double) * filter_order, "Butterworth A"));
return bw_coeff;
}
void ED_anim_free_butterworth_coefficients(ButterworthCoefficients *bw_coeff)
{
MEM_freeN(bw_coeff->d1);
MEM_freeN(bw_coeff->d2);
MEM_freeN(bw_coeff->A);
MEM_freeN(bw_coeff);
}
void ED_anim_calculate_butterworth_coefficients(const float cutoff_frequency,
const float sampling_frequency,
ButterworthCoefficients *bw_coeff)
{
double s = (double)sampling_frequency;
const double a = tan(M_PI * cutoff_frequency / s);
const double a2 = a * a;
double r;
for (int i = 0; i < bw_coeff->filter_order; ++i) {
r = sin(M_PI * (2.0 * i + 1.0) / (4.0 * bw_coeff->filter_order));
s = a2 + 2.0 * a * r + 1.0;
bw_coeff->A[i] = a2 / s;
bw_coeff->d1[i] = 2.0 * (1 - a2) / s;
bw_coeff->d2[i] = -(a2 - 2.0 * a * r + 1.0) / s;
}
}
static double butterworth_filter_value(
double x, double *w0, double *w1, double *w2, ButterworthCoefficients *bw_coeff)
{
for (int i = 0; i < bw_coeff->filter_order; i++) {
w0[i] = bw_coeff->d1[i] * w1[i] + bw_coeff->d2[i] * w2[i] + x;
x = bw_coeff->A[i] * (w0[i] + 2.0 * w1[i] + w2[i]);
w2[i] = w1[i];
w1[i] = w0[i];
}
return x;
}
static float butterworth_calculate_blend_value(float *samples,
float *filtered_values,
const int start_index,
const int end_index,
const int sample_index,
const int blend_in_out)
{
if (start_index == end_index || blend_in_out == 0) {
return samples[start_index];
}
const float blend_in_y_samples = samples[start_index];
const float blend_out_y_samples = samples[end_index];
const float blend_in_y_filtered = filtered_values[start_index + blend_in_out];
const float blend_out_y_filtered = filtered_values[end_index - blend_in_out];
const float slope_in_samples = samples[start_index] - samples[start_index - 1];
const float slope_out_samples = samples[end_index] - samples[end_index + 1];
const float slope_in_filtered = filtered_values[start_index + blend_in_out - 1] -
filtered_values[start_index + blend_in_out];
const float slope_out_filtered = filtered_values[end_index - blend_in_out] -
filtered_values[end_index - blend_in_out - 1];
if (sample_index - start_index <= blend_in_out) {
const int blend_index = sample_index - start_index;
const float blend_in_out_factor = clamp_f((float)blend_index / blend_in_out, 0.0f, 1.0f);
const float blend_value = interpf(blend_in_y_filtered +
slope_in_filtered * (blend_in_out - blend_index),
blend_in_y_samples + slope_in_samples * blend_index,
blend_in_out_factor);
return blend_value;
}
if (end_index - sample_index <= blend_in_out) {
const int blend_index = end_index - sample_index;
const float blend_in_out_factor = clamp_f((float)blend_index / blend_in_out, 0.0f, 1.0f);
const float blend_value = interpf(blend_out_y_filtered +
slope_out_filtered * (blend_in_out - blend_index),
blend_out_y_samples + slope_out_samples * blend_index,
blend_in_out_factor);
return blend_value;
}
return 0;
}
/**
* \param samples are expected to start at the first frame of the segment with a buffer of size
* segment->filter_order at the left.
*/
void butterworth_smooth_fcurve_segment(FCurve *fcu,
FCurveSegment *segment,
float *samples,
const int sample_count,
const float factor,
const int blend_in_out,
const int sample_rate,
ButterworthCoefficients *bw_coeff)
{
const int filter_order = bw_coeff->filter_order;
float *filtered_values = static_cast<float *>(
MEM_callocN(sizeof(float) * sample_count, "Butterworth Filtered FCurve Values"));
double *w0 = static_cast<double *>(MEM_callocN(sizeof(double) * filter_order, "w0"));
double *w1 = static_cast<double *>(MEM_callocN(sizeof(double) * filter_order, "w1"));
double *w2 = static_cast<double *>(MEM_callocN(sizeof(double) * filter_order, "w2"));
/* The values need to be offset so the first sample starts at 0. This avoids oscillations at the
* start and end of the curve. */
const float fwd_offset = samples[0];
for (int i = 0; i < sample_count; i++) {
const double x = (double)(samples[i] - fwd_offset);
const double filtered_value = butterworth_filter_value(x, w0, w1, w2, bw_coeff);
filtered_values[i] = (float)(filtered_value) + fwd_offset;
}
for (int i = 0; i < filter_order; i++) {
w0[i] = 0.0;
w1[i] = 0.0;
w2[i] = 0.0;
}
const float bwd_offset = filtered_values[sample_count - 1];
/* Run the filter backwards as well to remove phase offset. */
for (int i = sample_count - 1; i >= 0; i--) {
const double x = (double)(filtered_values[i] - bwd_offset);
const double filtered_value = butterworth_filter_value(x, w0, w1, w2, bw_coeff);
filtered_values[i] = (float)(filtered_value) + bwd_offset;
}
const int segment_end_index = segment->start_index + segment->length;
BezTriple left_bezt = fcu->bezt[segment->start_index];
BezTriple right_bezt = fcu->bezt[segment_end_index - 1];
const int samples_start_index = filter_order * sample_rate;
const int samples_end_index = (int)(right_bezt.vec[1][0] - left_bezt.vec[1][0] + filter_order) *
sample_rate;
const int blend_in_out_clamped = min_ii(blend_in_out,
(samples_end_index - samples_start_index) / 2);
for (int i = segment->start_index; i < segment_end_index; i++) {
float blend_in_out_factor;
if (blend_in_out_clamped == 0) {
blend_in_out_factor = 1;
}
else if (i < segment->start_index + segment->length / 2) {
blend_in_out_factor = min_ff((float)(i - segment->start_index) / blend_in_out_clamped, 1.0f);
}
else {
blend_in_out_factor = min_ff((float)(segment_end_index - i - 1) / blend_in_out_clamped,
1.0f);
}
const float x_delta = fcu->bezt[i].vec[1][0] - left_bezt.vec[1][0] + filter_order;
const int filter_index = (int)(x_delta * sample_rate);
const float blend_value = butterworth_calculate_blend_value(samples,
filtered_values,
samples_start_index,
samples_end_index,
filter_index,
blend_in_out_clamped);
const float blended_value = interpf(
filtered_values[filter_index], blend_value, blend_in_out_factor);
const float key_y_value = interpf(blended_value, samples[filter_index], factor);
BKE_fcurve_keyframe_move_value_with_handles(&fcu->bezt[i], key_y_value);
}
MEM_freeN(filtered_values);
MEM_freeN(w0);
MEM_freeN(w1);
MEM_freeN(w2);
}
/* ---------------- */
void ED_ANIM_get_1d_gauss_kernel(const float sigma, const int kernel_size, double *r_kernel)
@ -741,11 +932,13 @@ struct TempFrameValCache {
void sample_fcurve_segment(FCurve *fcu,
const float start_frame,
const int sample_rate,
float *samples,
const int sample_count)
{
for (int i = 0; i < sample_count; i++) {
samples[i] = evaluate_fcurve(fcu, start_frame + i);
const float evaluation_time = start_frame + ((float)i / sample_rate);
samples[i] = evaluate_fcurve(fcu, evaluation_time);
}
}

View File

@ -428,12 +428,32 @@ void blend_to_neighbor_fcurve_segment(struct FCurve *fcu,
struct FCurveSegment *segment,
float factor);
void breakdown_fcurve_segment(struct FCurve *fcu, struct FCurveSegment *segment, float factor);
/**
* Get a 1D gauss kernel. Since the kernel is symmetrical, only calculates the positive side.
* \param sigma: The shape of the gauss distribution.
* \param kernel_size: How long the kernel array is.
*/
void ED_ANIM_get_1d_gauss_kernel(const float sigma, int kernel_size, double *r_kernel);
typedef struct ButterworthCoefficients {
double *A, *d1, *d2;
int filter_order;
} ButterworthCoefficients;
ButterworthCoefficients *ED_anim_allocate_butterworth_coefficients(const int filter_order);
void ED_anim_free_butterworth_coefficients(struct ButterworthCoefficients *bw_coeff);
void ED_anim_calculate_butterworth_coefficients(float cutoff,
float sampling_frequency,
struct ButterworthCoefficients *bw_coeff);
void butterworth_smooth_fcurve_segment(struct FCurve *fcu,
struct FCurveSegment *segment,
float *samples,
int sample_count,
float factor,
int blend_in_out,
int sample_rate,
struct ButterworthCoefficients *bw_coeff);
void smooth_fcurve_segment(struct FCurve *fcu,
struct FCurveSegment *segment,
float *samples,
@ -452,10 +472,9 @@ void blend_to_default_fcurve(struct PointerRNA *id_ptr, struct FCurve *fcu, floa
*/
void smooth_fcurve(struct FCurve *fcu);
void sample_fcurve(struct FCurve *fcu);
void sample_fcurve_segment(struct FCurve *fcu,
float start_frame,
float *r_samples,
int sample_count);
/** \param sample_rate indicates how many samples per frame should be generated. */
void sample_fcurve_segment(
struct FCurve *fcu, float start_frame, int sample_rate, float *r_samples, int sample_count);
/* ----------- */

View File

@ -121,6 +121,7 @@ void GRAPH_OT_breakdown(struct wmOperatorType *ot);
void GRAPH_OT_ease(struct wmOperatorType *ot);
void GRAPH_OT_decimate(struct wmOperatorType *ot);
void GRAPH_OT_blend_to_default(struct wmOperatorType *ot);
void GRAPH_OT_butterworth_smooth(struct wmOperatorType *ot);
void GRAPH_OT_gaussian_smooth(struct wmOperatorType *ot);
void GRAPH_OT_sample(struct wmOperatorType *ot);
void GRAPH_OT_bake(struct wmOperatorType *ot);

View File

@ -473,6 +473,7 @@ void graphedit_operatortypes()
WM_operatortype_append(GRAPH_OT_ease);
WM_operatortype_append(GRAPH_OT_blend_to_default);
WM_operatortype_append(GRAPH_OT_gaussian_smooth);
WM_operatortype_append(GRAPH_OT_butterworth_smooth);
WM_operatortype_append(GRAPH_OT_euler_filter);
WM_operatortype_append(GRAPH_OT_delete);
WM_operatortype_append(GRAPH_OT_duplicate);

View File

@ -18,6 +18,7 @@
#include "MEM_guardedalloc.h"
#include "BLI_listbase.h"
#include "BLI_math.h"
#include "BLI_string.h"
#include "DNA_anim_types.h"
@ -1004,6 +1005,7 @@ struct tFCurveSegmentLink {
FCurve *fcu;
FCurveSegment *segment;
float *samples; /* Array of y-values of the FCurve segment. */
int sample_count;
};
static void gaussian_smooth_allocate_operator_data(tGraphSliderOp *gso,
@ -1037,7 +1039,7 @@ static void gaussian_smooth_allocate_operator_data(tGraphSliderOp *gso,
(filter_width * 2 + 1);
float *samples = static_cast<float *>(
MEM_callocN(sizeof(float) * sample_count, "Smooth FCurve Op Samples"));
sample_fcurve_segment(fcu, left_bezt.vec[1][0] - filter_width, samples, sample_count);
sample_fcurve_segment(fcu, left_bezt.vec[1][0] - filter_width, 1, samples, sample_count);
segment_link->samples = samples;
BLI_addtail(&segment_links, segment_link);
}
@ -1139,7 +1141,7 @@ static void gaussian_smooth_graph_keys(bAnimContext *ac,
(filter_width * 2 + 1);
float *samples = static_cast<float *>(
MEM_callocN(sizeof(float) * sample_count, "Smooth FCurve Op Samples"));
sample_fcurve_segment(fcu, left_bezt.vec[1][0] - filter_width, samples, sample_count);
sample_fcurve_segment(fcu, left_bezt.vec[1][0] - filter_width, 1, samples, sample_count);
smooth_fcurve_segment(fcu, segment, samples, factor, filter_width, kernel);
MEM_freeN(samples);
}
@ -1223,3 +1225,291 @@ void GRAPH_OT_gaussian_smooth(wmOperatorType *ot)
32);
}
/** \} */
/* -------------------------------------------------------------------- */
/** \name Butterworth Smooth Operator
* \{ */
typedef struct tBtwOperatorData {
ButterworthCoefficients *coefficients;
ListBase segment_links; /* tFCurveSegmentLink */
ListBase anim_data; /* bAnimListElem */
} tBtwOperatorData;
static int btw_calculate_sample_count(BezTriple *right_bezt,
BezTriple *left_bezt,
const int filter_order,
const int samples_per_frame)
{
/* Adding a constant 60 frames to combat the issue that the phase delay is shifting data out of
* the sample count range. This becomes an issue when running the filter backwards. */
const int sample_count = ((int)(right_bezt->vec[1][0] - left_bezt->vec[1][0]) + 1 +
(filter_order * 2)) *
samples_per_frame +
60;
return sample_count;
}
static void btw_smooth_allocate_operator_data(tGraphSliderOp *gso,
const int filter_order,
const int samples_per_frame)
{
tBtwOperatorData *operator_data = static_cast<tBtwOperatorData *>(
MEM_callocN(sizeof(tBtwOperatorData), "tBtwOperatorData"));
operator_data->coefficients = ED_anim_allocate_butterworth_coefficients(filter_order);
ListBase anim_data = {NULL, NULL};
ANIM_animdata_filter(
&gso->ac, &anim_data, OPERATOR_DATA_FILTER, gso->ac.data, eAnimCont_Types(gso->ac.datatype));
ListBase segment_links = {NULL, NULL};
LISTBASE_FOREACH (bAnimListElem *, ale, &anim_data) {
FCurve *fcu = (FCurve *)ale->key_data;
ListBase fcu_segments = find_fcurve_segments(fcu);
LISTBASE_FOREACH (FCurveSegment *, segment, &fcu_segments) {
tFCurveSegmentLink *segment_link = static_cast<tFCurveSegmentLink *>(
MEM_callocN(sizeof(tFCurveSegmentLink), "FCurve Segment Link"));
segment_link->fcu = fcu;
segment_link->segment = segment;
BezTriple left_bezt = fcu->bezt[segment->start_index];
BezTriple right_bezt = fcu->bezt[segment->start_index + segment->length - 1];
const int sample_count = btw_calculate_sample_count(
&right_bezt, &left_bezt, filter_order, samples_per_frame);
float *samples = static_cast<float *>(
MEM_callocN(sizeof(float) * sample_count, "Btw Smooth FCurve Op Samples"));
sample_fcurve_segment(
fcu, left_bezt.vec[1][0] - filter_order, samples_per_frame, samples, sample_count);
segment_link->samples = samples;
segment_link->sample_count = sample_count;
BLI_addtail(&segment_links, segment_link);
}
}
operator_data->anim_data = anim_data;
operator_data->segment_links = segment_links;
gso->operator_data = operator_data;
}
static void btw_smooth_free_operator_data(void *operator_data)
{
tBtwOperatorData *btw_data = (tBtwOperatorData *)operator_data;
LISTBASE_FOREACH (tFCurveSegmentLink *, segment_link, &btw_data->segment_links) {
MEM_freeN(segment_link->samples);
MEM_freeN(segment_link->segment);
}
ED_anim_free_butterworth_coefficients(btw_data->coefficients);
BLI_freelistN(&btw_data->segment_links);
ANIM_animdata_freelist(&btw_data->anim_data);
MEM_freeN(btw_data);
}
static void btw_smooth_modal_update(bContext *C, wmOperator *op)
{
tGraphSliderOp *gso = static_cast<tGraphSliderOp *>(op->customdata);
bAnimContext ac;
if (ANIM_animdata_get_context(C, &ac) == 0) {
return;
}
common_draw_status_header(C, gso, "Butterworth Smooth");
tBtwOperatorData *operator_data = (tBtwOperatorData *)gso->operator_data;
const float frame_rate = (float)(ac.scene->r.frs_sec) / ac.scene->r.frs_sec_base;
const int samples_per_frame = RNA_int_get(op->ptr, "samples_per_frame");
const float sampling_frequency = frame_rate * samples_per_frame;
const float cutoff_frequency = slider_factor_get_and_remember(op);
const int blend_in_out = RNA_int_get(op->ptr, "blend_in_out");
ED_anim_calculate_butterworth_coefficients(
cutoff_frequency, sampling_frequency, operator_data->coefficients);
LISTBASE_FOREACH (tFCurveSegmentLink *, segment, &operator_data->segment_links) {
butterworth_smooth_fcurve_segment(segment->fcu,
segment->segment,
segment->samples,
segment->sample_count,
1,
blend_in_out,
samples_per_frame,
operator_data->coefficients);
}
LISTBASE_FOREACH (bAnimListElem *, ale, &operator_data->anim_data) {
ale->update |= ANIM_UPDATE_DEFAULT;
}
ANIM_animdata_update(&ac, &operator_data->anim_data);
WM_event_add_notifier(C, NC_ANIMATION | ND_KEYFRAME | NA_EDITED, NULL);
}
static int btw_smooth_invoke(bContext *C, wmOperator *op, const wmEvent *event)
{
const int invoke_result = graph_slider_invoke(C, op, event);
if (invoke_result == OPERATOR_CANCELLED) {
return invoke_result;
}
tGraphSliderOp *gso = static_cast<tGraphSliderOp *>(op->customdata);
gso->modal_update = btw_smooth_modal_update;
gso->factor_prop = RNA_struct_find_property(op->ptr, "cutoff_frequency");
const int filter_order = RNA_int_get(op->ptr, "filter_order");
const int samples_per_frame = RNA_int_get(op->ptr, "samples_per_frame");
btw_smooth_allocate_operator_data(gso, filter_order, samples_per_frame);
gso->free_operator_data = btw_smooth_free_operator_data;
const float frame_rate = (float)(gso->scene->r.frs_sec) / gso->scene->r.frs_sec_base;
const float sampling_frequency = frame_rate * samples_per_frame;
ED_slider_factor_bounds_set(gso->slider, 0, sampling_frequency / 2);
ED_slider_factor_set(gso->slider, RNA_float_get(op->ptr, "cutoff_frequency"));
ED_slider_allow_overshoot_set(gso->slider, false, false);
ED_slider_mode_set(gso->slider, SLIDER_MODE_FLOAT);
ED_slider_unit_set(gso->slider, "Hz");
common_draw_status_header(C, gso, "Butterworth Smooth");
return invoke_result;
}
static void btw_smooth_graph_keys(bAnimContext *ac,
const float factor,
const int blend_in_out,
float cutoff_frequency,
const int filter_order,
const int samples_per_frame)
{
ListBase anim_data = {NULL, NULL};
ANIM_animdata_filter(
ac, &anim_data, OPERATOR_DATA_FILTER, ac->data, eAnimCont_Types(ac->datatype));
ButterworthCoefficients *bw_coeff = ED_anim_allocate_butterworth_coefficients(filter_order);
const float frame_rate = (float)(ac->scene->r.frs_sec) / ac->scene->r.frs_sec_base;
const float sampling_frequency = frame_rate * samples_per_frame;
/* Clamp cutoff frequency to Nyquist Frequency. */
cutoff_frequency = min_ff(cutoff_frequency, sampling_frequency / 2);
ED_anim_calculate_butterworth_coefficients(cutoff_frequency, sampling_frequency, bw_coeff);
LISTBASE_FOREACH (bAnimListElem *, ale, &anim_data) {
FCurve *fcu = (FCurve *)ale->key_data;
ListBase segments = find_fcurve_segments(fcu);
LISTBASE_FOREACH (FCurveSegment *, segment, &segments) {
BezTriple left_bezt = fcu->bezt[segment->start_index];
BezTriple right_bezt = fcu->bezt[segment->start_index + segment->length - 1];
const int sample_count = btw_calculate_sample_count(
&right_bezt, &left_bezt, filter_order, samples_per_frame);
float *samples = static_cast<float *>(
MEM_callocN(sizeof(float) * sample_count, "Smooth FCurve Op Samples"));
sample_fcurve_segment(
fcu, left_bezt.vec[1][0] - filter_order, samples_per_frame, samples, sample_count);
butterworth_smooth_fcurve_segment(
fcu, segment, samples, sample_count, factor, blend_in_out, samples_per_frame, bw_coeff);
MEM_freeN(samples);
}
BLI_freelistN(&segments);
ale->update |= ANIM_UPDATE_DEFAULT;
}
ED_anim_free_butterworth_coefficients(bw_coeff);
ANIM_animdata_update(ac, &anim_data);
ANIM_animdata_freelist(&anim_data);
}
static int btw_smooth_exec(bContext *C, wmOperator *op)
{
bAnimContext ac;
if (ANIM_animdata_get_context(C, &ac) == 0) {
return OPERATOR_CANCELLED;
}
const float blend = RNA_float_get(op->ptr, "blend");
const float cutoff_frequency = RNA_float_get(op->ptr, "cutoff_frequency");
const int filter_order = RNA_int_get(op->ptr, "filter_order");
const int samples_per_frame = RNA_int_get(op->ptr, "samples_per_frame");
const int blend_in_out = RNA_int_get(op->ptr, "blend_in_out");
btw_smooth_graph_keys(
&ac, blend, blend_in_out, cutoff_frequency, filter_order, samples_per_frame);
/* Set notifier that keyframes have changed. */
WM_event_add_notifier(C, NC_ANIMATION | ND_KEYFRAME | NA_EDITED, NULL);
return OPERATOR_FINISHED;
}
void GRAPH_OT_butterworth_smooth(wmOperatorType *ot)
{
/* Identifiers. */
ot->name = "Butterworth Smooth";
ot->idname = "GRAPH_OT_butterworth_smooth";
ot->description = "Smooth an F-Curve while maintaining the general shape of the curve";
/* API callbacks. */
ot->invoke = btw_smooth_invoke;
ot->modal = graph_slider_modal;
ot->exec = btw_smooth_exec;
ot->poll = graphop_editable_keyframes_poll;
/* Flags. */
ot->flag = OPTYPE_REGISTER | OPTYPE_UNDO | OPTYPE_BLOCKING | OPTYPE_GRAB_CURSOR_X;
RNA_def_float(ot->srna,
"cutoff_frequency",
3.0f,
0.0f,
FLT_MAX,
"Frquency Cutoff (Hz)",
"Lower values give a smoother curve",
0.0f,
FLT_MAX);
RNA_def_int(ot->srna,
"filter_order",
4,
1,
32,
"Filter Order",
"Higher values produce a harder frequency cutoff",
1,
16);
RNA_def_int(ot->srna,
"samples_per_frame",
1,
1,
64,
"Samples per Frame",
"How many samples to calculate per frame, helps with subframe data",
1,
16);
RNA_def_float_factor(ot->srna,
"blend",
1.0f,
0,
FLT_MAX,
"Blend",
"How much to blend to the smoothed curve",
0.0f,
1.0f);
RNA_def_int(ot->srna,
"blend_in_out",
1,
0,
INT_MAX,
"Blend In/Out",
"Linearly blend the smooth data to the border frames of the selection",
0,
128);
}
/** \} */