diff --git a/tools/.gitignore b/tools/.gitignore
index ac00e6bcd2..3073e2e521 100644
--- a/tools/.gitignore
+++ b/tools/.gitignore
@@ -1 +1,2 @@
 xpubaddrgen
+mksecptable
diff --git a/tools/Makefile b/tools/Makefile
index 67c652637f..519bb1de9a 100644
--- a/tools/Makefile
+++ b/tools/Makefile
@@ -23,7 +23,7 @@ CFLAGS   += $(OPTFLAGS) \
             -Werror \
             -I..
 
-all: xpubaddrgen
+all: xpubaddrgen mksecptable
 
 OBJS = ../bip32.o ../ecdsa.o ../sha2.o ../bignum.o ../base58.o ../secp256k1.o ../ripemd160.o ../hmac.o ../rand.o
 
@@ -33,5 +33,8 @@ OBJS = ../bip32.o ../ecdsa.o ../sha2.o ../bignum.o ../base58.o ../secp256k1.o ..
 xpubaddrgen: xpubaddrgen.o $(OBJS)
 	$(CC) xpubaddrgen.o $(OBJS) -o xpubaddrgen
 
+mksecptable: mksecptable.o $(OBJS)
+	$(CC) mksecptable.o $(OBJS) -o mksecptable
+
 clean:
-	rm -f *.o xpubaddrgen
+	rm -f *.o xpubaddrgen mksecptable
diff --git a/tools/README.md b/tools/README.md
index 8f24f19ab0..93ea9fb1de 100644
--- a/tools/README.md
+++ b/tools/README.md
@@ -43,3 +43,11 @@ Example output:
 It will print ```<jobid> error``` when there was an error processing job jobid.
 
 It will print ```error``` when it encountered a malformed line.
+
+
+mksecptable
+-----------
+
+mksecptable computes the points of the form `(2*j+1)*16^i*G` and prints them in the format to be included in `secp256k1.c`.  These points are used by the fast ECC multiplication.
+
+It is only meant to be run if the `scalar_mult` algorithm changes.
diff --git a/tools/mksecptable.c b/tools/mksecptable.c
new file mode 100644
index 0000000000..e6a1997d90
--- /dev/null
+++ b/tools/mksecptable.c
@@ -0,0 +1,59 @@
+#include <stdio.h>
+#include <assert.h>
+#include "bignum.h"
+#include "ecdsa.h"
+#include "secp256k1.h"
+#include "rand.h"
+
+/*
+ * This program prints the contents of the secp256k1_cp array.
+ * The entry secp256k1_cp[i][j] contains the number (2*j+1)*16^i*G,
+ * where G is the generator of secp256k1.
+ */
+int main(int __attribute__((unused)) argc, char __attribute__((unused)) **argv) {
+	int i,j,k;
+	curve_point ng = G256k1;
+	curve_point pow2ig = G256k1;
+#ifndef NDEBUG
+	init_rand(); // needed for point_multiply()
+#endif
+	for (i = 0; i < 64; i++) {
+		// invariants:
+		//   pow2ig = 16^i * G
+		//   ng     = pow2ig
+		printf("\t{\n");
+		for (j = 0; j < 8; j++) {
+			// invariants:
+			//   pow2ig = 16^i * G
+			//   ng     = (2*j+1) * 16^i * G
+#ifndef NDEBUG
+			curve_point checkresult;
+			bignum256 a;
+			bn_zero(&a);
+			a.val[(4*i) / 30] = ((uint32_t) 2*j+1) << ((4*i) % 30);
+			bn_normalize(&a);
+			point_multiply(&a, &G256k1, &checkresult);
+			assert(point_is_equal(&checkresult, &ng));
+#endif
+			printf("\t\t/* %2d*16^%d*G: */\n\t\t{{{", 2*j + 1, i);
+			// print x coordinate
+			for (k = 0; k < 9; k++) {
+				printf((k < 8 ? "0x%08x, " : "0x%04x"), ng.x.val[k]);
+			}
+			printf("}},\n\t\t {{");
+			// print y coordinate
+			for (k = 0; k < 9; k++) {
+				printf((k < 8 ? "0x%08x, " : "0x%04x"), ng.y.val[k]);
+			}
+			if (j == 7) {
+				printf("}}}\n\t},\n");
+			} else {
+				printf("}}},\n");
+				point_add(&pow2ig, &ng);
+			}
+			point_add(&pow2ig, &ng);
+		}
+		pow2ig = ng;
+	}
+	return 0;
+}