bmp.c

#include <stdio.h>

#include "bmp.h"

typedef struct __attribute__((packed)) {
	uint16_t bfType;  // specifies the file type
	uint32_t bfSize;  // specifies the size in bytes of the bitmap file
	uint16_t bfReserved1;  // reserved; must be 0
	uint16_t bfReserved2;  // reserved; must be 0
	uint32_t bfOffBits;  // species the offset in bytes from the bitmapfileheader to the bitmap bits
} BMPFileHeader;

void *LoadGrayscaleBitmapFile(const char *filename, BMPInfoHeader *bitmapInfoHeader)
{
	// open filename in read binary mode
	FILE *filePtr = fopen(filename,"rb");
	if (filePtr == NULL) {
		return NULL;
	}

	// read the bitmap file header
	BMPFileHeader bitmapFileHeader;
	if (fread(&bitmapFileHeader, sizeof bitmapFileHeader, 1, filePtr) != 1) {
		fclose(filePtr);
		return NULL;
	}

	// verify that this is a bmp file by checking the bitmap ID
	if (bitmapFileHeader.bfType != 0x4D42) {
		fclose(filePtr);
		return NULL;
	}

	// read the bitmap info header
	if (fread(bitmapInfoHeader, sizeof *bitmapInfoHeader, 1, filePtr) != 1) {
		fclose(filePtr);
		return NULL;
	}

	// move cursor to the beginning of bitmap data
	fseek(filePtr, bitmapFileHeader.bfOffBits, SEEK_SET);

	// allocate enough memory for the bitmap image data
	uint8_t *bitmapImage = malloc(bitmapInfoHeader->biSizeImage);

	// read in the bitmap image data
	if (fread(bitmapImage, bitmapInfoHeader->biSizeImage, 1, filePtr) != 1) {
		fclose(filePtr);
		free(bitmapImage);
		return NULL;
	}

	double *tmp = malloc(bitmapInfoHeader->biSizeImage / 4);

	for (uint32_t i = 0; i < bitmapInfoHeader->biSizeImage / 4; i++) {
		tmp[i] = bitmapImage[i * 4] / 255.0;
	}

	// close file and return bitmap iamge data
	fclose(filePtr);

	free(bitmapImage);
	return tmp;
}

bmp.h

#ifndef BMP_H
#define BMP_H

#include <stdint.h>
#include <stdlib.h>

typedef struct __attribute__((packed))
{
	uint32_t biSize;  // specifies the number of bytes required by the struct
	int32_t  biWidth;  // specifies width in pixels
	int32_t  biHeight;  // species height in pixels
	uint16_t biPlanes; // specifies the number of color planes, must be 1
	uint16_t biBitCount; // specifies the number of bit per pixel
	uint32_t biCompression;// spcifies the type of compression
	uint32_t biSizeImage;  // size of image in bytes
	int32_t  biXPelsPerMeter;  // number of pixels per meter in x axis
	int32_t  biYPelsPerMeter;  // number of pixels per meter in y axis
	uint32_t biClrUsed;  // number of colors used by th ebitmap
	uint32_t biClrImportant;  // number of colors that are important
} BMPInfoHeader;

void *LoadGrayscaleBitmapFile(const char *filename, BMPInfoHeader *bitmapInfoHeader);

#endif // BMP_H

week2.c

#include <stdio.h>
#include <assert.h>
#include <math.h>

#include "bmp.h"

void filter3(
	uint32_t width,
	uint32_t height,
	const double (*indata)[width],
	double (**outdata)[width]
)
{
	double (*tmp)[width] = malloc(height * sizeof tmp[0]);

	for (int32_t r = 0; r < height; r++) {
		for (int32_t c = 0; c < width; c++) {

			double coeff = 1.0 / 9;

			// "replicate" boundary condition
			int32_t rm = r - 1 < 0 ? 0 : r - 1;
			int32_t rp = r + 1 >= height ? height - 1 : r + 1;
			int32_t cm = c - 1 < 0 ? 0 : c - 1;
			int32_t cp = c + 1 >= width ? width - 1 : c + 1;

			tmp[r][c] = coeff * (
				indata[rm][cm] +
				indata[rm][c ] +
				indata[rm][cp] +
				indata[r ][cm] +
				indata[r ][c ] +
				indata[r ][cp] +
				indata[rp][cm] +
				indata[rp][c ] +
				indata[rp][cp]
			);
		}
	}

	*outdata = tmp;
}

void filter5(
	uint32_t width,
	uint32_t height,
	const double (*indata)[width],
	double (**outdata)[width]
)
{
	double (*tmp)[width] = malloc(height * sizeof tmp[0]);

	for (int32_t r = 0; r < height; r++) {
		for (int32_t c = 0; c < width; c++) {

			double coeff = 1.0 / 25;

			// "replicate" boundary condition
			int32_t rm2 = r - 2 < 0 ? 0 : r - 2;
			int32_t rm1 = r - 1 < 0 ? 0 : r - 1;
			int32_t rp1 = r + 1 >= height ? height - 1 : r + 1;
			int32_t rp2 = r + 2 >= height ? height - 1 : r + 2;

			int32_t cm2 = c - 2 < 0 ? 0 : c - 2;
			int32_t cm1 = c - 1 < 0 ? 0 : c - 1;
			int32_t cp1 = c + 1 >= width ? width - 1 : c + 1;
			int32_t cp2 = c + 2 >= width ? width - 1 : c + 2;

			// What's a loop???
			tmp[r][c] = coeff * (
				indata[rm2][cm2] +
				indata[rm2][cm1] +
				indata[rm2][c  ] +
				indata[rm2][cp1] +
				indata[rm2][cp2] +

				indata[rm1][cm2] +
				indata[rm1][cm1] +
				indata[rm1][c  ] +
				indata[rm1][cp1] +
				indata[rm1][cp2] +

				indata[r  ][cm2] +
				indata[r  ][cm1] +
				indata[r  ][c  ] +
				indata[r  ][cp1] +
				indata[r  ][cp2] +

				indata[rp1][cm2] +
				indata[rp1][cm1] +
				indata[rp1][c  ] +
				indata[rp1][cp1] +
				indata[rp1][cp2] +

				indata[rp2][cm2] +
				indata[rp2][cm1] +
				indata[rp2][c  ] +
				indata[rp2][cp1] +
				indata[rp2][cp2]
			);
		}
	}

	*outdata = tmp;
}

double PSNR(
	uint32_t width,
	uint32_t height,
	const double (*x)[width],
	const double (*y)[width]
)
{
	double maxI = 1.0;
	double mse = 0.0;

	for (uint32_t r = 0; r < height; r++) {
		for (uint32_t c = 0; c < width; c++) {
			mse += (x[r][c] - y[r][c]) * (x[r][c] - y[r][c]);
		}
	}

	mse *= 1.0 / (width * height);

	return 10 * log10(maxI * maxI / mse);
}

int main()
{
	BMPInfoHeader hdr;
	void *tmp = LoadGrayscaleBitmapFile("week2_quiz.bmp", &hdr);
	double (*in)[hdr.biWidth] = tmp;

	double (*out)[hdr.biWidth] = NULL;
	filter3(hdr.biWidth, hdr.biHeight, in, &out); // or filter5

	double psnr = PSNR(hdr.biWidth, hdr.biHeight, in, out);
	printf("%lf\n", psnr);

	free(in);
	free(out);
}