diff --git a/README.md b/README.md
index 24e2fc5..6b4a9e2 100644
--- a/README.md
+++ b/README.md
@@ -28,6 +28,18 @@ Implemented are the following [finite impulse response](https://en.wikipedia.org
 * [high-pass filter](https://en.wikipedia.org/wiki/High-pass_filter)
 * [band-pass filter](https://en.wikipedia.org/wiki/Band-pass_filter)
 
+### [ema.hh](ema.hh)
+
+The [exponential moving average](https://en.wikipedia.org/wiki/Moving_average#Exponential_moving_average) is an [infinite impulse response](https://en.wikipedia.org/wiki/Infinite_impulse_response) [low-pass filter](https://en.wikipedia.org/wiki/Low-pass_filter).
+There is also support for cascading, to improve [roll-off](https://en.wikipedia.org/wiki/Roll-off) while a correction factor helps to keep the same [cutoff frequency](https://en.wikipedia.org/wiki/Cutoff_frequency).
+
+### [biquad.hh](biquad.hh)
+
+The following [infinite impulse response](https://en.wikipedia.org/wiki/Infinite_impulse_response) [digital biquad filter](https://en.wikipedia.org/wiki/Digital_biquad_filter) implementations are available:
+
+* [second-order Butterworth low pass filter](https://en.wikipedia.org/wiki/Butterworth_filter)
+* [2n-order Butterworth cascade of second-order low pass filters](https://en.wikipedia.org/wiki/Butterworth_filter)
+
 ### [const.hh](const.hh)
 
 Some constants we need
diff --git a/biquad.hh b/biquad.hh
new file mode 100644
index 0000000..e2db34a
--- /dev/null
+++ b/biquad.hh
@@ -0,0 +1,71 @@
+/*
+Digital biquad filter
+
+Copyright 2019 Ahmet Inan <inan@aicodix.de>
+*/
+
+#pragma once
+
+#include "const.hh"
+#include "unit_circle.hh"
+
+namespace DSP {
+
+template <typename TYPE, typename VALUE>
+class Biquad
+{
+	TYPE x1, x2, y1, y2;
+	VALUE b0a0, b1a0, b2a0, a1a0, a2a0;
+public:
+	constexpr Biquad() : x1(0), x2(0), y1(0), y2(0),
+		b0a0(1), b1a0(0), b2a0(0), a1a0(0), a2a0(0)
+	{
+	}
+	void lowpass(int n, int N, VALUE Q = Const<VALUE>::InvSqrtTwo())
+	{
+		VALUE alpha = UnitCircle<VALUE>::sin(n, N) / (VALUE(2) * Q),
+			cn = UnitCircle<VALUE>::cos(n, N),
+			b0 = (VALUE(1) - cn) / VALUE(2),
+			b1 = VALUE(1) - cn,
+			b2 = (VALUE(1) - cn) / VALUE(2),
+			a0 = VALUE(1) + alpha,
+			a1 = -VALUE(2) * cn,
+			a2 = VALUE(1) - alpha;
+		b0a0 = b0 / a0;
+		b1a0 = b1 / a0;
+		b2a0 = b2 / a0;
+		a1a0 = a1 / a0;
+		a2a0 = a2 / a0;
+	}
+	TYPE operator()(TYPE x0)
+	{
+		TYPE y0 = b0a0*x0 + b1a0*x1 + b2a0*x2
+			- a1a0*y1 - a2a0*y2;
+		x2 = x1; x1 = x0;
+		y2 = y1; y1 = y0;
+		return y0;
+	}
+};
+
+template <typename TYPE, typename VALUE, int NUM = 1>
+class BiquadCascade
+{
+	static const int ORDER = 2 * NUM;
+	Biquad<TYPE, VALUE> cascade[NUM];
+public:
+	void lowpass(int n, int N)
+	{
+		for (int i = 0; i < NUM; ++i)
+			cascade[i].lowpass(n, N, VALUE(1) / (VALUE(2) *
+				UnitCircle<VALUE>::cos(2*i+1, 4*ORDER)));
+	}
+	TYPE operator()(TYPE input)
+	{
+		for (int i = 0; i < NUM; ++i)
+			input = cascade[i](input);
+		return input;
+	}
+};
+
+}
+
diff --git a/ema.hh b/ema.hh
new file mode 100644
index 0000000..210d699
--- /dev/null
+++ b/ema.hh
@@ -0,0 +1,88 @@
+/*
+Exponential moving average
+
+Copyright 2019 Ahmet Inan <inan@aicodix.de>
+*/
+
+#include "utils.hh"
+#include "unit_circle.hh"
+
+#pragma once
+
+namespace DSP {
+
+template <typename TYPE, typename VALUE>
+class EMA
+{
+	TYPE prev;
+	VALUE alpha;
+	static VALUE cutoff_alpha(int n, int N)
+	{
+		VALUE x = UnitCircle<VALUE>::cos(n, N);
+		return x-VALUE(1)+sqrt(x*(x-VALUE(4))+VALUE(3));
+	}
+public:
+	EMA() : prev(0), alpha(1)
+	{
+	}
+	void cutoff(int n, int N)
+	{
+		alpha = cutoff_alpha(n, N);
+	}
+	TYPE operator()(TYPE input)
+	{
+		return prev = lerp(alpha, prev, input);
+	}
+};
+
+template <typename TYPE, typename VALUE, int ORDER = 1>
+class EMACascade
+{
+	TYPE prev[ORDER];
+	VALUE alpha;
+	static constexpr VALUE cascade_factor(int n)
+	{
+		// sqrt(e(1/n*l(2))-1)
+		return
+			n == 2 ? .6435942529055826247354434374182098089242027424440076511561520093520L :
+			n == 3 ? .5098245285339585980866134620415975171004662971729213347866667412748L :
+			n == 4 ? .4349794420460822959023617403105260887838801269492020945122292703957L :
+			n == 5 ? .3856142567346739149847213965799930737314794827929705449459805476567L :
+			n == 6 ? .3499457790992384314823879950215014587750945565870355960474164484675L :
+			n == 7 ? .3226290651410878970708971817580301448745259798477682170908870012270L :
+			n == 8 ? .3008450309798346315492310951214492771360484456002599192820777432131L :
+			n == 9 ? .2829482972069387977799123584440613415063420296903386876512703264517L :
+			n == 10 ? .2679056970956257290860485974661891235831324413057223464282157983007L :
+			n == 11 ? .2550315459702243937103098780342370981449477556708628034249027758407L :
+			n == 12 ? .2438505574307659124524083486415529630644763241448445300317390218311L :
+			n == 13 ? .2340215299532217691188557628487431016851279721568792383213019527102L :
+			n == 14 ? .2252923404228812762935182605055552181633447787610798690912140601327L :
+			n == 15 ? .2174721196397982365782144895603276387766971846333067133969338153904L :
+			n == 16 ? .2104133608576552353576644943161359426249259240013118976757629219108L :
+			1;
+	}
+	static VALUE cutoff_alpha(int n, int N)
+	{
+		VALUE x = UnitCircle<VALUE>::cos(n, N);
+		return x-VALUE(1)+sqrt(x*(x-VALUE(4))+VALUE(3));
+	}
+public:
+	EMACascade() : alpha(1)
+	{
+		for (int i = 0; i < ORDER; ++i)
+			prev[i] = TYPE(0);
+	}
+	void cutoff(int n, int N)
+	{
+		alpha = cutoff_alpha(nearbyint(n / cascade_factor(ORDER)), N);
+	}
+	TYPE operator()(TYPE input)
+	{
+		for (int i = 0; i < ORDER; ++i)
+			prev[i] = input = lerp(alpha, prev[i], input);
+		return input;
+	}
+};
+
+}
+