diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml
index 1c3ba83..cbf2da2 100644
--- a/.github/workflows/ci.yml
+++ b/.github/workflows/ci.yml
@@ -27,7 +27,7 @@ jobs:
     needs: [get_notebooks]
     strategy:
       matrix:
-        python-version: [3.8, 3.9, '3.10', '3.11', '3.12']
+        python-version: [3.9, '3.10', '3.11', '3.12', '3.13', '3.14']
         notebook: ${{fromJson(needs.get_notebooks.outputs.notebook)}}
       fail-fast: false
     name: Execute notebooks
@@ -53,7 +53,7 @@ jobs:
     strategy:
       matrix:
         os: [ubuntu-latest, macos-latest, windows-latest]
-        python-version: [3.8, 3.9, '3.10', '3.11', '3.12']
+        python-version: [3.9, '3.10', '3.11', '3.12', '3.13', '3.14']
       fail-fast: false
     name: Run tests
     steps:
diff --git a/README.md b/README.md
index 0c0816a..03047ad 100644
--- a/README.md
+++ b/README.md
@@ -3,7 +3,7 @@
 [![PyPI Version](https://img.shields.io/pypi/v/pyrtools.svg)](https://pypi.org/project/pyrtools/)
 [![Anaconda-Server Badge](https://anaconda.org/conda-forge/pyrtools/badges/version.svg)](https://anaconda.org/conda-forge/pyrtools)
 [![License: MIT](https://img.shields.io/badge/License-MIT-yellow.svg)](https://github.com/LabForComputationalVision/pyrtools/blob/main/LICENSE)
-![Python version](https://img.shields.io/badge/python-3.8|3.9|3.10|3.11|3.12-blue.svg)
+![Python version](https://img.shields.io/badge/python-3.9|3.10|3.11|3.12|3.13|3.14-blue.svg)
 [![Build Status](https://github.com/LabForComputationalVision/pyrtools/workflows/build/badge.svg)](https://github.com/LabForComputationalVision/pyrtools/actions?query=workflow%3Abuild)
 [![Documentation Status](https://readthedocs.org/projects/pyrtools/badge/?version=latest)](https://pyrtools.readthedocs.io/en/latest/?badge=latest)
 [![DOI](https://zenodo.org/badge/137527035.svg)](https://zenodo.org/doi/10.5281/zenodo.10161031)
diff --git a/docs/index.rst b/docs/index.rst
index 84d5871..d05587b 100644
--- a/docs/index.rst
+++ b/docs/index.rst
@@ -7,7 +7,7 @@
 .. |license-shield| image:: https://img.shields.io/badge/license-MIT-yellow.svg
 			    :target: https://github.com/LabForComputationalVision/pyrtools/blob/main/LICENSE
 
-.. |python-version-shield| image:: https://img.shields.io/badge/python-3.8%7C3.9%7C3.10%7C3.11%7C3.12-blue.svg
+.. |python-version-shield| image:: https://img.shields.io/badge/python-3.9%7C3.10%7C3.11%7C3.12%7C3.13%7C3.14-blue.svg
 
 .. |build| image:: https://github.com/LabForComputationalVision/pyrtools/workflows/build/badge.svg
 		     :target: https://github.com/LabForComputationalVision/pyrtools/actions?query=workflow%3Abuild
diff --git a/pyproject.toml b/pyproject.toml
index 18f9d67..d469b17 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -4,11 +4,11 @@ dynamic = ["version"]
 authors = [{name="Pyrtools authors"}]
 description = "Python tools for multi-scale image processing, including Laplacian pyramids, Wavelets, and Steerable Pyramids."
 readme = "README.md"
-requires-python = ">=3.8"
+license = "MIT"
+requires-python = ">=3.9"
 classifiers = [
     "Development Status :: 4 - Beta",
     "Programming Language :: Python :: 3",
-    "License :: OSI Approved :: MIT License",
     "Intended Audience :: Science/Research",
 ]
 keywords = ['image processing', 'visual information processing', 'computational models']
@@ -34,7 +34,7 @@ docs = [
 ]
 
 [build-system]
-requires = ["setuptools", "wheel", "setuptools-scm[toml]"]
+requires = ["setuptools>=70.1", "setuptools-scm[toml]"]
 build-backend = "setuptools.build_meta"
 
 [project.urls]
diff --git a/setup.py b/setup.py
index 34545e0..5a197cd 100644
--- a/setup.py
+++ b/setup.py
@@ -1,7 +1,7 @@
 #! /usr/bin/env python
 
-from wheel.bdist_wheel import bdist_wheel
 from setuptools import setup, Extension
+from setuptools.command.bdist_wheel import bdist_wheel
 
 # Adapted from the cibuildwheel example https://github.com/joerick/python-ctypes-package-sample
 # it marks the wheel as not specific to the Python API version.
diff --git a/src/pyrtools/pyramids/c/convolve.c b/src/pyrtools/pyramids/c/convolve.c
index 3d4e94d..1341a6f 100755
--- a/src/pyrtools/pyramids/c/convolve.c
+++ b/src/pyrtools/pyramids/c/convolve.c
@@ -69,19 +69,10 @@
 	}
 */
 
-int internal_reduce(image, x_dim, y_dim, filt, temp, x_fdim, y_fdim,
-		x_start, x_step, x_stop, y_start, y_step, y_stop,
-		result, edges)
-  register image_type *image, *temp;
-  register int x_fdim, x_dim;
-  register image_type *result;
-  register int x_step, y_step;
-  int x_start, y_start;
-  int x_stop, y_stop;     
-  image_type *filt; 
-  int y_dim, y_fdim;
-  char *edges;
-  { 
+int internal_reduce(register image_type *image, register int x_dim, int y_dim, image_type *filt, register image_type *temp, register int x_fdim, int y_fdim,
+		int x_start, register int x_step, int x_stop, int y_start, register int y_step, int y_stop,
+		register image_type *result, char *edges)
+  {
   register double sum;
   register int filt_pos, im_pos, x_filt_stop;
   register int x_pos, filt_size = x_fdim*y_fdim;
@@ -236,18 +227,9 @@ int internal_reduce(image, x_dim, y_dim, filt, temp, x_fdim, y_fdim,
       }									\
   }
 
-int internal_expand(image,filt,temp,x_fdim,y_fdim,
-		x_start,x_step,x_stop,y_start,y_step,y_stop,
-		result,x_dim,y_dim,edges)
-  register image_type *result, *temp;
-  register int x_fdim, x_dim;
-  register int x_step, y_step;
-  register image_type *image; 
-  int x_start, y_start;
-  int x_stop, y_stop;
-  image_type *filt; 
-  int y_fdim, y_dim;
-  char *edges;
+int internal_expand(register image_type *image, image_type *filt, register image_type *temp, register int x_fdim, int y_fdim,
+		int x_start, register int x_step, int x_stop, int y_start, register int y_step, int y_stop,
+		register image_type *result, register int x_dim, int y_dim, char *edges)
   {
   register double val;
   register int filt_pos, res_pos, x_filt_stop;
diff --git a/src/pyrtools/pyramids/c/convolve.h b/src/pyrtools/pyramids/c/convolve.h
index 9f15987..27f619c 100755
--- a/src/pyrtools/pyramids/c/convolve.h
+++ b/src/pyrtools/pyramids/c/convolve.h
@@ -24,7 +24,7 @@
 #define AND &&
 #define OR ||
 
-typedef  int (*fptr)();
+typedef int (*fptr)(double *, int , int , int , int , double *, int );
 
 typedef struct 
   {
diff --git a/src/pyrtools/pyramids/c/edges.c b/src/pyrtools/pyramids/c/edges.c
index 1ecdf8d..d322fbe 100755
--- a/src/pyrtools/pyramids/c/edges.c
+++ b/src/pyrtools/pyramids/c/edges.c
@@ -40,8 +40,15 @@ below.
 #define sgn(a)  ( ((a)>0)?1:(((a)<0)?-1:0) )
 #define clip(a,mn,mx)  ( ((a)<(mn))?(mn):(((a)>=(mx))?(mx-1):(a)) )
 
-int reflect1(), reflect2(), qreflect2(), repeat(), zero(), Extend(), nocompute();
-int ereflect(), predict();
+int nocompute(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
+int zero(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
+int reflect1(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
+int reflect2(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
+int qreflect2(register double *filt, int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
+int repeat(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
+int Extend(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
+int predict(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
+int ereflect(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e);
 
 /* Lookup table matching a descriptive string to the edge-handling function */
 #if !THINK_C
@@ -138,10 +145,7 @@ r_or_e - equal to one of the two constants EXPAND or REDUCE.
 nocompute() - Return zero for values where filter hangs over the edge.
 */
 
-int nocompute(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  register double *filt, *result;
-  register int x_dim;
-  int y_dim, x_pos, y_pos, r_or_e;
+int nocompute(register double *filt, register int x_dim, int y_dim, int x_pos, int y_pos, register double *result, int r_or_e)
   {
   register int i;
   register int size = x_dim*y_dim;
@@ -157,10 +161,7 @@ int nocompute(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
 /* --------------------------------------------------------------------
 zero() - Zero outside of image.  Discontinuous, but adds zero energy. */
 
-int zero(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  register double *filt, *result;
-  register int x_dim;
-  int y_dim, x_pos, y_pos, r_or_e;
+int zero(double *filt, int x_dim, int y_dim, int x_pos, int y_pos, double *result, int r_or_e)
   {
   register int y_filt,x_filt, y_res,x_res;
   int filt_sz = x_dim*y_dim;
@@ -190,10 +191,7 @@ are subsampling by 2, since it maintains parity (even pixels positions
 remain even, odd ones remain odd).
 */	 
 
-int reflect1(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  register double *filt, *result;
-  register int x_dim;
-  int y_dim, x_pos, y_pos, r_or_e;
+int reflect1(double *filt, int x_dim, int y_dim, int x_pos, int y_pos, double *result, int r_or_e)
   {
   int filt_sz = x_dim*y_dim;
   register int y_filt,x_filt, y_res, x_res;
@@ -259,10 +257,7 @@ then the next pixel, etc.  Continuous, but discontinuous first
 derivative.
 */
 
-int reflect2(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  register double *filt, *result;
-  register int x_dim;
-  int y_dim, x_pos, y_pos, r_or_e;
+int reflect2(double *filt, int x_dim, int y_dim, int x_pos, int y_pos, double *result, int r_or_e)
   {
   int filt_sz = x_dim*y_dim;
   register int y_filt,x_filt, y_res, x_res;
@@ -334,9 +329,7 @@ qreflect2() - Modified version of reflect2 that  works properly for
 even-length QMF filters.
 */
 
-int qreflect2(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  double *filt, *result;
-  int x_dim, y_dim, x_pos, y_pos, r_or_e;
+int qreflect2(double *filt, int x_dim, int y_dim, int x_pos, int y_pos, double *result, int r_or_e)
   {
   reflect2(filt,x_dim,y_dim,x_pos,y_pos,result,0);
   return(0);
@@ -347,10 +340,7 @@ repeat() - repeat edge pixel.  Continuous, with discontinuous first
 derivative.
 */
 
-int repeat(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  register double *filt, *result;
-  register int x_dim;
-  int y_dim, x_pos, y_pos, r_or_e;
+int repeat(double *filt, int x_dim, int y_dim, int x_pos, int y_pos, double *result, int r_or_e)
   {
   register int y_filt,x_filt, y_res,x_res, y_tmp, x_tmp;
   register int x_base = (x_pos>0)?(x_dim-1):0;
@@ -418,10 +408,7 @@ value.  Maintains continuity in intensity AND first derivative (but
 not higher derivs).
 */
 
-int Extend(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  register double *filt, *result;
-  register int x_dim;
-  int y_dim, x_pos, y_pos, r_or_e;
+int Extend(double *filt, int x_dim, int y_dim, int x_pos, int y_pos, double *result, int r_or_e)
   {
   int filt_sz = x_dim*y_dim;
   register int y_filt,x_filt, y_res,x_res, y_tmp, x_tmp;
@@ -538,10 +525,7 @@ by the reciprocal of the percentage of filter being used.  (i.e. if
 50% of the filter is hanging over the edge of the image, multiply the
 taps being used by 2).  */
 
-int predict(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  register double *filt, *result;
-  register int x_dim;
-  int y_dim, x_pos, y_pos, r_or_e;
+int predict(double *filt, int x_dim, int y_dim, int x_pos, int y_pos, double *result, int r_or_e)
   {
   register int y_filt,x_filt, y_res,x_res;
   register double taps_used = 0.0; /* int *** */
@@ -583,10 +567,7 @@ by root 2.  This maintains orthogonality of odd-length linear-phase
 QMF filters, but it is not useful for most applications, since it
 alters the DC level.  */
 
-int ereflect(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
-  register double *filt, *result;
-  register int x_dim;
-  int y_dim, x_pos, y_pos, r_or_e;
+int ereflect(double *filt, int x_dim, int y_dim, int x_pos, int y_pos, double *result, int r_or_e)
   {
   register int y_filt,x_filt, y_res,x_res;
   register int x_base = (x_pos>0)?(x_dim-1):0;
diff --git a/src/pyrtools/pyramids/c/internal_pointOp.c b/src/pyrtools/pyramids/c/internal_pointOp.c
index 314fd81..5579752 100644
--- a/src/pyrtools/pyramids/c/internal_pointOp.c
+++ b/src/pyrtools/pyramids/c/internal_pointOp.c
@@ -5,10 +5,7 @@
 /* Use linear interpolation on a lookup table.
    Taken from OBVIUS.  EPS, Spring, 1987.
  */
-void internal_pointop (im, res, size, lut, lutsize, origin, increment, warnings)
-  register double *im, *res, *lut;
-  register double origin, increment; 
-  register int size, lutsize, warnings;
+void internal_pointop (register double *im, register double *res, register int size, register double *lut, register int lutsize, register double origin, register double increment, register int warnings)
   {
   register int i, index;
   register double pos;
diff --git a/src/pyrtools/pyramids/c/wrap.c b/src/pyrtools/pyramids/c/wrap.c
index 3651b6b..4599199 100755
--- a/src/pyrtools/pyramids/c/wrap.c
+++ b/src/pyrtools/pyramids/c/wrap.c
@@ -173,13 +173,9 @@ int internal_wrap_reduce(image, x_dim, y_dim, filt, x_fdim, y_fdim,
 	  imval[YIND][XIND] += val * filt[filt_pos]; \
       }
 
-int internal_wrap_expand(image, filt, x_fdim, y_fdim,
-	      x_start, x_step, x_stop, y_start, y_step, y_stop,
-	      result, x_dim, y_dim)
-  register image_type *filt, *result;
-  register int x_fdim, y_fdim, x_dim, y_dim;
-  image_type *image; 
-  int x_start, x_step, x_stop, y_start, y_step, y_stop;
+int internal_wrap_expand(image_type *image, register image_type *filt, register int x_fdim, register int y_fdim,
+	      int x_start, int x_step, int x_stop, int y_start, int y_step, int y_stop,
+	      register image_type *result, register int x_dim, register int y_dim)
   {
   register double val;
   register int filt_size = x_fdim*y_fdim;