Two-sided exponential impulse response blur

#define EXP_BLUR_ROW_PASSES 2
#define EXP_BLUR_COL_PASSES 2
// Returns +1 for even integers and -1 for odd integers
#define EXP_BLUR_DIRECTION(X) ((1 ^ (X)) - (X))

static inline void blur_inner(
	unsigned char *pixel,
	int pixel_size,
	int *zR,
	int *zG,
	int *zB,
	int *zA,
	int alpha,
	int aprec,
	int zprec)
{
	*zR += (alpha * ((*(pixel    ) << zprec) - *zR)) >> aprec;
	*zG += (alpha * ((*(pixel + 1) << zprec) - *zG)) >> aprec;
	*zB += (alpha * ((*(pixel + 2) << zprec) - *zB)) >> aprec;

	*(pixel    ) = *zR >> zprec;
	*(pixel + 1) = *zG >> zprec;
	*(pixel + 2) = *zB >> zprec;

	if (pixel_size > 3) {
		*zA += (alpha * ((*(pixel + 3) << zprec) - *zA)) >> aprec;
		*(pixel + 3) = *zA >> zprec;
	}
}

static inline void blur_row(
	unsigned char *pixels,
	int width,
	int height,
	int pixel_size,
	int row_size,
	int y,
	int alpha,
	int aprec,
	int zprec)
{
	unsigned char *row = &pixels[y * row_size];

	int zR = *(row    ) << zprec;
	int zG = *(row + 1) << zprec;
	int zB = *(row + 2) << zprec;
	int zA = pixel_size > 3 ? *(row + 3) << zprec : 0;

	for (int pass = 0; pass < EXP_BLUR_ROW_PASSES; ++pass)
		for (int index = 0; index < width; ++index)
			blur_inner(
				&row[(index * EXP_BLUR_DIRECTION(pass) + (pass % 2) * (width - 1)) * pixel_size],
				pixel_size,
				&zR,
				&zG,
				&zB,
				&zA,
				alpha,
				aprec,
				zprec
			);
}

static inline void blur_col(
	unsigned char *pixels,
	int width,
	int height,
	int pixel_size,
	int row_size,
	int x,
	int alpha,
	int aprec,
	int zprec)
{
	unsigned char *col = &pixels[x * pixel_size];

	int zR = *(col    ) << zprec;
	int zG = *(col + 1) << zprec;
	int zB = *(col + 2) << zprec;
	int zA = pixel_size > 3 ? *(col + 3) << zprec : 0;

	for (int pass = 0; pass < EXP_BLUR_COL_PASSES; ++pass)
		for (int index = 0; index < height; ++index)
			blur_inner(
				&col[(index * EXP_BLUR_DIRECTION(pass) + (pass % 2) * (height - 1)) * row_size],
				pixel_size,
				&zR,
				&zG,
				&zB,
				&zA,
				alpha,
				aprec,
				zprec
			);
}

static void exp_blur(
	unsigned char *pixels, // image data
	int width,             // image width in pixels
	int height,            // image height in pixels
	int pixel_size,        // pixel size in bytes
	double radius,         // kernel radius
	int aprec,             // precision of alpha parameter in fixed-point format 0.aprec
	int zprec)             // precision of state parameters zR,zG,zB and zA in fp format 8.zprec
{
	int row_size = width * pixel_size;
	if (aprec < 0) {
		aprec = 16;
	}
	if (zprec < 0) {
		zprec = 7;
	}

	// Calculate the alpha such that 90% of the kernel is within the radius (kernel extends to infinity)
	int alpha = (int) ((1 << aprec) * (1.0f - expf(-2.3f / (radius + 1.0f))));

	for (int row = 0; row < height; ++row)
		blur_row(
			pixels,
			width,
			height,
			pixel_size,
			row_size,
			row,
			alpha,
			aprec,
			zprec
		);

	for (int col = 0; col < width; ++col)
		blur_col(
			pixels,
			width,
			height,
			pixel_size,
			row_size,
			col,
			alpha,
			aprec,
			zprec
		);
}