diff --git a/decode.cc b/decode.cc
index 9a1b8a1..8d601f3 100644
--- a/decode.cc
+++ b/decode.cc
@@ -28,6 +28,7 @@ namespace DSP { using std::abs; using std::min; using std::cos; using std::sin;
 #include "crc.hh"
 #include "osd.hh"
 #include "psk.hh"
+#include "qam.hh"
 #include "polar_tables.hh"
 #include "polar_helper.hh"
 #include "polar_encoder.hh"
@@ -172,7 +173,7 @@ struct Decoder
 #endif
 	typedef DSP::Const<value> Const;
 	static const int code_order = 11;
-	static const int mod_bits = 2;
+	static const int mod_bits = 4;
 	static const int code_len = 1 << code_order;
 	static const int symbol_len = (1280 * rate) / 8000;
 	static const int filter_len = (((21 * rate) / 8000) & ~3) | 1;
@@ -180,7 +181,7 @@ struct Decoder
 	static const int extended_len = symbol_len + guard_len;
 	static const int max_bits = 1360 + 32;
 	static const int cons_cols = 256;
-	static const int cons_rows = 4;
+	static const int cons_rows = 2;
 	static const int cons_total = cons_rows * cons_cols;
 	static const int code_off = - cons_cols / 2;
 	static const int mls0_len = 127;
@@ -252,15 +253,15 @@ struct Decoder
 	}
 	cmplx mod_map(code_type *b)
 	{
-		return PhaseShiftKeying<4, cmplx, code_type>::map(b);
+		return QuadratureAmplitudeModulation<16, cmplx, code_type>::map(b);
 	}
 	void mod_hard(code_type *b, cmplx c)
 	{
-		PhaseShiftKeying<4, cmplx, code_type>::hard(b, c);
+		QuadratureAmplitudeModulation<16, cmplx, code_type>::hard(b, c);
 	}
 	void mod_soft(code_type *b, cmplx c, value precision)
 	{
-		PhaseShiftKeying<4, cmplx, code_type>::soft(b, c, precision);
+		QuadratureAmplitudeModulation<16, cmplx, code_type>::soft(b, c, precision);
 	}
 	const cmplx *next_sample()
 	{
@@ -415,7 +416,7 @@ struct Decoder
 				value snr = DSP::decibel(precision);
 				std::cerr << " " << snr;
 				for (int i = 0; i < cons_cols; ++i)
-					mod_soft(code+2*(cons_cols*j+i), cons[cons_cols*j+i], precision);
+					mod_soft(code+mod_bits*(cons_cols*j+i), cons[cons_cols*j+i], precision);
 			}
 			std::cerr << std::endl;
 		} else {
diff --git a/encode.cc b/encode.cc
index 6277b58..7d947e2 100644
--- a/encode.cc
+++ b/encode.cc
@@ -18,6 +18,7 @@ Copyright 2021 Ahmet Inan <inan@aicodix.de>
 #include "mls.hh"
 #include "crc.hh"
 #include "psk.hh"
+#include "qam.hh"
 #include "polar_tables.hh"
 #include "polar_helper.hh"
 #include "polar_encoder.hh"
@@ -27,14 +28,14 @@ template <typename value, typename cmplx, int rate>
 struct Encoder
 {
 	typedef int8_t code_type;
-	static const int mod_bits = 2;
+	static const int mod_bits = 4;
 	static const int code_order = 11;
 	static const int code_len = 1 << code_order;
 	static const int symbol_len = (1280 * rate) / 8000;
 	static const int guard_len = symbol_len / 8;
 	static const int max_bits = 1360 + 32;
 	static const int cons_cols = 256;
-	static const int cons_rows = 4;
+	static const int cons_rows = 2;
 	static const int mls0_len = 127;
 	static const int mls0_poly = 0b10001001;
 	static const int mls1_len = 255;
@@ -187,7 +188,7 @@ struct Encoder
 	}
 	cmplx mod_map(code_type *b)
 	{
-		return PhaseShiftKeying<4, cmplx, code_type>::map(b);
+		return QuadratureAmplitudeModulation<16, cmplx, code_type>::map(b);
 	}
 	Encoder(DSP::WritePCM<value> *pcm, const uint8_t *inp, int freq_off, uint64_t call_sign, int oper_mode, int reserved_tones) :
 		pcm(pcm), crc0(0xA8F4), crc1(0x8F6E37A0), bchenc({
diff --git a/qam.hh b/qam.hh
new file mode 100644
index 0000000..7715e3e
--- /dev/null
+++ b/qam.hh
@@ -0,0 +1,254 @@
+/*
+Quadrature amplitude modulation
+
+Copyright 2018 Ahmet Inan <xdsopl@gmail.com>
+*/
+
+#ifndef QAM_HH
+#define QAM_HH
+
+template <int NUM, typename TYPE, typename CODE>
+struct QuadratureAmplitudeModulation;
+
+template <typename TYPE, typename CODE>
+struct QuadratureAmplitudeModulation<16, TYPE, CODE>
+{
+	static const int NUM = 16;
+	static const int BITS = 4;
+	typedef TYPE complex_type;
+	typedef typename TYPE::value_type value_type;
+	typedef CODE code_type;
+
+	static constexpr value_type FAC = 1.0540925533894596;
+	static constexpr value_type RCP = 3 * FAC;
+	static constexpr value_type AMP = 1 / RCP;
+	static constexpr value_type DIST = 2 * AMP;
+
+	static constexpr value_type amp(int i)
+	{
+		return AMP * i;
+	}
+
+	static code_type quantize(value_type precision, value_type value)
+	{
+		value *= DIST * precision;
+		if (std::is_integral<code_type>::value)
+			value = std::nearbyint(value);
+		if (std::is_same<code_type, int8_t>::value)
+			value = std::min<value_type>(std::max<value_type>(value, -127), 127);
+		return value;
+	}
+
+	static void hard(code_type *b, complex_type c)
+	{
+		b[0] = c.real() < amp(0) ? code_type(-1) : code_type(1);
+		b[1] = c.imag() < amp(0) ? code_type(-1) : code_type(1);
+		b[2] = std::abs(c.real()) < amp(2) ? code_type(-1) : code_type(1);
+		b[3] = std::abs(c.imag()) < amp(2) ? code_type(-1) : code_type(1);
+	}
+
+	static void soft(code_type *b, complex_type c, value_type precision)
+	{
+		b[0] = quantize(precision, c.real());
+		b[1] = quantize(precision, c.imag());
+		b[2] = quantize(precision, std::abs(c.real())-amp(2));
+		b[3] = quantize(precision, std::abs(c.imag())-amp(2));
+	}
+
+	static complex_type map(code_type *b)
+	{
+		return AMP * complex_type(
+			b[0]*(b[2]+value_type(2)),
+			b[1]*(b[3]+value_type(2))
+		);
+	}
+};
+
+template <typename TYPE, typename CODE>
+struct QuadratureAmplitudeModulation<64, TYPE, CODE>
+{
+	static const int NUM = 64;
+	static const int BITS = 6;
+	typedef TYPE complex_type;
+	typedef typename TYPE::value_type value_type;
+	typedef CODE code_type;
+
+	static constexpr value_type FAC = 0.9258200997725516;
+	static constexpr value_type RCP = 7 * FAC;
+	static constexpr value_type AMP = 1 / RCP;
+	static constexpr value_type DIST = 2 * AMP;
+
+	static constexpr value_type amp(int i)
+	{
+		return AMP * i;
+	}
+
+	static code_type quantize(value_type precision, value_type value)
+	{
+		value *= DIST * precision;
+		if (std::is_integral<code_type>::value)
+			value = std::nearbyint(value);
+		if (std::is_same<code_type, int8_t>::value)
+			value = std::min<value_type>(std::max<value_type>(value, -127), 127);
+		return value;
+	}
+
+	static void hard(code_type *b, complex_type c)
+	{
+		b[0] = c.real() < amp(0) ? code_type(-1) : code_type(1);
+		b[1] = c.imag() < amp(0) ? code_type(-1) : code_type(1);
+		b[2] = std::abs(c.real()) < amp(4) ? code_type(-1) : code_type(1);
+		b[3] = std::abs(c.imag()) < amp(4) ? code_type(-1) : code_type(1);
+		b[4] = std::abs(std::abs(c.real())-amp(4)) < amp(2) ? code_type(-1) : code_type(1);
+		b[5] = std::abs(std::abs(c.imag())-amp(4)) < amp(2) ? code_type(-1) : code_type(1);
+	}
+
+	static void soft(code_type *b, complex_type c, value_type precision)
+	{
+		b[0] = quantize(precision, c.real());
+		b[1] = quantize(precision, c.imag());
+		b[2] = quantize(precision, std::abs(c.real())-amp(4));
+		b[3] = quantize(precision, std::abs(c.imag())-amp(4));
+		b[4] = quantize(precision, std::abs(std::abs(c.real())-amp(4))-amp(2));
+		b[5] = quantize(precision, std::abs(std::abs(c.imag())-amp(4))-amp(2));
+	}
+
+	static complex_type map(code_type *b)
+	{
+		return AMP * complex_type(
+			b[0]*(b[2]*(b[4]+value_type(2))+value_type(4)),
+			b[1]*(b[3]*(b[5]+value_type(2))+value_type(4))
+		);
+	}
+};
+
+template <typename TYPE, typename CODE>
+struct QuadratureAmplitudeModulation<256, TYPE, CODE>
+{
+	static const int NUM = 256;
+	static const int BITS = 8;
+	typedef TYPE complex_type;
+	typedef typename TYPE::value_type value_type;
+	typedef CODE code_type;
+
+	static constexpr value_type FAC = 0.8692269873603529;
+	static constexpr value_type RCP = 15 * FAC;
+	static constexpr value_type AMP = 1 / RCP;
+	static constexpr value_type DIST = 2 * AMP;
+
+	static constexpr value_type amp(int i)
+	{
+		return AMP * i;
+	}
+
+	static code_type quantize(value_type precision, value_type value)
+	{
+		value *= DIST * precision;
+		if (std::is_integral<code_type>::value)
+			value = std::nearbyint(value);
+		if (std::is_same<code_type, int8_t>::value)
+			value = std::min<value_type>(std::max<value_type>(value, -127), 127);
+		return value;
+	}
+
+	static void hard(code_type *b, complex_type c)
+	{
+		b[0] = c.real() < amp(0) ? code_type(-1) : code_type(1);
+		b[1] = c.imag() < amp(0) ? code_type(-1) : code_type(1);
+		b[2] = std::abs(c.real()) < amp(8) ? code_type(-1) : code_type(1);
+		b[3] = std::abs(c.imag()) < amp(8) ? code_type(-1) : code_type(1);
+		b[4] = std::abs(std::abs(c.real())-amp(8)) < amp(4) ? code_type(-1) : code_type(1);
+		b[5] = std::abs(std::abs(c.imag())-amp(8)) < amp(4) ? code_type(-1) : code_type(1);
+		b[6] = std::abs(std::abs(std::abs(c.real())-amp(8))-amp(4)) < amp(2) ? code_type(-1) : code_type(1);
+		b[7] = std::abs(std::abs(std::abs(c.imag())-amp(8))-amp(4)) < amp(2) ? code_type(-1) : code_type(1);
+	}
+
+	static void soft(code_type *b, complex_type c, value_type precision)
+	{
+		b[0] = quantize(precision, c.real());
+		b[1] = quantize(precision, c.imag());
+		b[2] = quantize(precision, std::abs(c.real())-amp(8));
+		b[3] = quantize(precision, std::abs(c.imag())-amp(8));
+		b[4] = quantize(precision, std::abs(std::abs(c.real())-amp(8))-amp(4));
+		b[5] = quantize(precision, std::abs(std::abs(c.imag())-amp(8))-amp(4));
+		b[6] = quantize(precision, std::abs(std::abs(std::abs(c.real())-amp(8))-amp(4))-amp(2));
+		b[7] = quantize(precision, std::abs(std::abs(std::abs(c.imag())-amp(8))-amp(4))-amp(2));
+	}
+
+	static complex_type map(code_type *b)
+	{
+		return AMP * complex_type(
+			b[0]*(b[2]*(b[4]*(b[6]+value_type(2))+value_type(4))+value_type(8)),
+			b[1]*(b[3]*(b[5]*(b[7]+value_type(2))+value_type(4))+value_type(8))
+		);
+	}
+};
+
+template <typename TYPE, typename CODE>
+struct QuadratureAmplitudeModulation<1024, TYPE, CODE>
+{
+	static const int NUM = 1024;
+	static const int BITS = 10;
+	typedef TYPE complex_type;
+	typedef typename TYPE::value_type value_type;
+	typedef CODE code_type;
+
+	static constexpr value_type FAC = 0.8424235391742344;
+	static constexpr value_type RCP = 31 * FAC;
+	static constexpr value_type AMP = 1 / RCP;
+	static constexpr value_type DIST = 2 * AMP;
+
+	static constexpr value_type amp(int i)
+	{
+		return AMP * i;
+	}
+
+	static code_type quantize(value_type precision, value_type value)
+	{
+		value *= DIST * precision;
+		if (std::is_integral<code_type>::value)
+			value = std::nearbyint(value);
+		if (std::is_same<code_type, int8_t>::value)
+			value = std::min<value_type>(std::max<value_type>(value, -127), 127);
+		return value;
+	}
+
+	static void hard(code_type *b, complex_type c)
+	{
+		b[0] = c.real() < amp(0) ? code_type(-1) : code_type(1);
+		b[1] = c.imag() < amp(0) ? code_type(-1) : code_type(1);
+		b[2] = std::abs(c.real()) < amp(16) ? code_type(-1) : code_type(1);
+		b[3] = std::abs(c.imag()) < amp(16) ? code_type(-1) : code_type(1);
+		b[4] = std::abs(std::abs(c.real())-amp(16)) < amp(8) ? code_type(-1) : code_type(1);
+		b[5] = std::abs(std::abs(c.imag())-amp(16)) < amp(8) ? code_type(-1) : code_type(1);
+		b[6] = std::abs(std::abs(std::abs(c.real())-amp(16))-amp(8)) < amp(4) ? code_type(-1) : code_type(1);
+		b[7] = std::abs(std::abs(std::abs(c.imag())-amp(16))-amp(8)) < amp(4) ? code_type(-1) : code_type(1);
+		b[8] = std::abs(std::abs(std::abs(std::abs(c.real())-amp(16))-amp(8))-amp(4)) < amp(2) ? code_type(-1) : code_type(1);
+		b[9] = std::abs(std::abs(std::abs(std::abs(c.imag())-amp(16))-amp(8))-amp(4)) < amp(2) ? code_type(-1) : code_type(1);
+	}
+
+	static void soft(code_type *b, complex_type c, value_type precision)
+	{
+		b[0] = quantize(precision, c.real());
+		b[1] = quantize(precision, c.imag());
+		b[2] = quantize(precision, std::abs(c.real())-amp(16));
+		b[3] = quantize(precision, std::abs(c.imag())-amp(16));
+		b[4] = quantize(precision, std::abs(std::abs(c.real())-amp(16))-amp(8));
+		b[5] = quantize(precision, std::abs(std::abs(c.imag())-amp(16))-amp(8));
+		b[6] = quantize(precision, std::abs(std::abs(std::abs(c.real())-amp(16))-amp(8))-amp(4));
+		b[7] = quantize(precision, std::abs(std::abs(std::abs(c.imag())-amp(16))-amp(8))-amp(4));
+		b[8] = quantize(precision, std::abs(std::abs(std::abs(std::abs(c.real())-amp(16))-amp(8))-amp(4))-amp(2));
+		b[9] = quantize(precision, std::abs(std::abs(std::abs(std::abs(c.imag())-amp(16))-amp(8))-amp(4))-amp(2));
+	}
+
+	static complex_type map(code_type *b)
+	{
+		return AMP * complex_type(
+			b[0]*(b[2]*(b[4]*(b[6]*(b[8]+value_type(2))+value_type(4))+value_type(8))+value_type(16)),
+			b[1]*(b[3]*(b[5]*(b[7]*(b[9]+value_type(2))+value_type(4))+value_type(8))+value_type(16))
+		);
+	}
+};
+
+#endif
+