#include #include #include #include #define IMAGE_WIDTH 43200 #define IMAGE_HEIGHT 16800 #define MASK_WIDTH 10800 #define MASK_HEIGHT 5400 #define OUTPUT_FILE_TEMPLATE "./Baseline_ETo_Data-Pass_%d.bin" #define FILENAME_MAX_LENGTH 40 #define HEADER_SIZE 32 long unsigned CROPPED_TOP_PIXELS = (MASK_WIDTH * MASK_HEIGHT * 10 / 180); // These will be set by findPixelRange(). uint16_t minPixelValue = 0; uint16_t maxPixelValue = 0xFFFF; double bitReductionFactor = 256; /** Copies the big-endian byte representation of the specified value into the specified buffer. */ void copyBytes(void* input, uint8_t* output, int unsigned length) { int unsigned isBigEndian = 1; isBigEndian = *((uint8_t*)(&isBigEndian)) == 0; for (int unsigned i = 0; i < length; i++) { int unsigned index = isBigEndian ? i : length - i - 1; output[i] = *((uint8_t*) input + index); } } /** * Write file header to the specified buffer. The header format is documented in the README. */ void setHeader(uint8_t *header) { for (int unsigned i = 0; i < HEADER_SIZE; i++) { header[i] = 0; } uint32_t width = IMAGE_WIDTH; uint32_t height = IMAGE_HEIGHT; // originally 0.1, then multiplied by a value to compensate for the bit depth reduction and divided by 25.4 to convert to inches. float scalingFactor = 0.1 * bitReductionFactor / 25.4; float minimumETo = minPixelValue * 0.1 / 25.4; // Version header[0] = 1; // Width copyBytes(&width, &(header[1]), 4); // Height copyBytes(&height, &(header[5]), 4); // Bit depth header[9] = 8; // Minimum ETo copyBytes(&minimumETo, &(header[10]), 4); // Scaling factor copyBytes(&scalingFactor, &(header[14]), 4); } /** * Calculates the minimum and maximum pixel values used in the image. These values can be used to optimally reduce the * bit depth by mapping the minimum value to 0 and the maximum value to 254 (reserving 255 for fill pixels) and linearly * interpolating the rest of the values. */ void findPixelRange(uint16_t* minPtr, uint16_t* maxPtr, double* bitReductionFactorPtr) { time_t startTime = clock(); uint16_t minValue = 0xFFFF; uint16_t maxValue = 0; FILE* inputFile = fopen("./MOD16A3_PET_2000_to_2013_mean.bin", "rb"); if (inputFile == NULL) { printf("An error occurred opening image file while finding min/max value.\n"); exit(1); } uint16_t buffer[IMAGE_WIDTH]; for (int unsigned y = 0; y < IMAGE_HEIGHT; y++) { if (y % 1000 == 0) { printf("Finding pixel range on row %d...\n", y); } fread(buffer, 2, IMAGE_WIDTH, inputFile); if (ferror(inputFile)) { printf("An error occurred reading image row %d while finding min/max values.\n", y); exit(1); } if (feof(inputFile)) { printf("Encountered EOF reading image row %d while finding min/max values.\n", y); exit(1); } for (unsigned int x = 0; x < IMAGE_WIDTH; x++) { uint16_t pixel = buffer[x]; // Skip fill pixels. if (pixel > 65528) { continue; } minValue = pixel < minValue ? pixel : minValue; maxValue = pixel > maxValue ? pixel : maxValue; } } *minPtr = minValue; *maxPtr = maxValue; *bitReductionFactorPtr = (maxValue - minValue + 1) / (float) 256; fclose(inputFile); printf("Found pixel range in %.1f seconds. Min value: %d\t Max value: %d\t Bit reduction factor:%f.\n", (clock() - startTime) / (float) CLOCKS_PER_SEC, minValue, maxValue, *bitReductionFactorPtr); } /** Reduces the image bit depth from 16 bits to 8 bits. */ void reduceBitDepth() { clock_t startTime = clock(); FILE* originalFile = fopen("./MOD16A3_PET_2000_to_2013_mean.bin", "rb"); if (originalFile == NULL) { printf("An error occurred opening input image file while reducing bit depth.\n"); exit(1); } char* reducedFileName = malloc(FILENAME_MAX_LENGTH); snprintf(reducedFileName, FILENAME_MAX_LENGTH, OUTPUT_FILE_TEMPLATE, 0); FILE* reducedFile = fopen(reducedFileName, "wb"); if (reducedFile == NULL) { printf("An error occurred opening output image file while reducing bit depth.\n"); exit(1); } // Write the file header. uint8_t header[32]; setHeader(header); fwrite(header, 1, 32, reducedFile); if (ferror(reducedFile)) { printf("An error occurred writing file header while reducing bit depth.\n"); exit(1); } uint16_t inputBuffer[IMAGE_WIDTH]; uint8_t outputBuffer[IMAGE_WIDTH]; for (int unsigned y = 0; y < IMAGE_HEIGHT; y++) { if (y % 1000 == 0) { printf("Reducing bit depth on row %d...\n", y); } fread(inputBuffer, 2, IMAGE_WIDTH, originalFile); if (ferror(originalFile)) { printf("An error occurred reading row %d while reducing bit depth.\n", y); exit(1); } if (feof(originalFile)) { printf("Encountered EOF reading row %d while reducing bit depth.\n", y); exit(1); } for (unsigned int x = 0; x < IMAGE_WIDTH; x++) { uint16_t originalPixel = inputBuffer[x]; uint8_t reducedPixel = originalPixel > 65528 ? 255 : (uint8_t) ((originalPixel - minPixelValue) / bitReductionFactor); outputBuffer[x] = reducedPixel; } fwrite(outputBuffer, 1, IMAGE_WIDTH, reducedFile); if (ferror(reducedFile)) { printf("An error occurred writing row %d while reducing bit depth.\n", y); exit(1); } } fclose(reducedFile); fclose(originalFile); printf("Finished reducing bit depth in %.1f seconds.\n", (clock() - startTime) / (double) CLOCKS_PER_SEC); } void fillMissingPixels(int unsigned pass) { clock_t startTime = clock(); char* inputFileName = malloc(FILENAME_MAX_LENGTH); snprintf(inputFileName, FILENAME_MAX_LENGTH, OUTPUT_FILE_TEMPLATE, pass - 1); FILE* inputFile = fopen(inputFileName, "rb"); if (inputFile == NULL) { printf("An error occurred opening input image file on pass %d.\n", pass); exit(1); } char* outputFileName = malloc(FILENAME_MAX_LENGTH); snprintf(outputFileName, FILENAME_MAX_LENGTH, OUTPUT_FILE_TEMPLATE, pass); FILE* outputFile = fopen(outputFileName, "wb"); if (outputFile == NULL) { printf("An error occurred opening output image file on pass %d.\n", pass); exit(1); } FILE* maskFile = fopen("./Ocean_Mask.bin", "rb"); if (maskFile == NULL) { printf("An error occurred opening mask image on pass %d.\n", pass); exit(1); } uint8_t outputBuffer[IMAGE_WIDTH]; // Skip the header. fseek(inputFile, 32, SEEK_SET); if (ferror(inputFile)) { printf("An error occurred reading header on pass %d.\n", pass); exit(1); } if (feof(inputFile)) { printf("Encountered EOF reading header on pass %d.\n", pass); exit(1); } // Write the file header. uint8_t header[32]; setHeader(header); fwrite(header, 1, 32, outputFile); if (ferror(outputFile)) { printf("An error occurred writing file header on pass %d.\n", pass); exit(1); } uint8_t* rows[5] = {0, 0, 0, 0, 0}; // Read the first 2 rows. for (int unsigned rowIndex = 3; rowIndex < 5; rowIndex++) { uint8_t* row = (uint8_t*) malloc(IMAGE_WIDTH); fread(row, 1, IMAGE_WIDTH, inputFile); if (ferror(inputFile)) { printf("An error occurred reading image row %d on pass %d.\n", rowIndex - 3, pass); exit(1); } if (feof(inputFile)) { printf("Encountered EOF reading image row %d on pass %d.\n", rowIndex - 3, pass); exit(1); } rows[rowIndex] = row; } long unsigned fixedPixels = 0; long unsigned unfixablePixels = 0; long unsigned waterPixels = 0; for (int unsigned y = 0; y < IMAGE_HEIGHT; y++) { if (y % 1000 == 0) { printf("Filling missing pixels on pass %d row %d.\n", pass, y); } // Read a row from the mask. uint8_t maskRow[MASK_WIDTH]; int unsigned maskOffset = y / (IMAGE_WIDTH / MASK_WIDTH) * MASK_WIDTH + CROPPED_TOP_PIXELS; fseek(maskFile, maskOffset, SEEK_SET); fread(maskRow, 1, MASK_WIDTH, maskFile); if (ferror(maskFile)) { printf("An error occurred reading mask at offset %d on pass %d.\n", maskOffset, pass); exit(1); } if (feof(maskFile)) { printf("Encountered EOF reading mask at offset %d on pass %d.\n", maskOffset, pass); exit(1); } // Free the oldest row. free(rows[0]); // Shift the previous rows back. for (int unsigned rowIndex = 1; rowIndex < 5; rowIndex++) { rows[rowIndex - 1] = rows[rowIndex]; } // Read the next row if one exists. if (y < IMAGE_HEIGHT - 2) { uint8_t* row = malloc(IMAGE_WIDTH); fread(row, 1, IMAGE_WIDTH, inputFile); if (ferror(inputFile)) { printf("An error occurred reading image row %d on pass %d.\n", y + 2, pass); exit(1); } if (feof(inputFile)) { printf("Encountered EOF reading image row %d on pass %d,\n", y + 2, pass); exit(1); } rows[4] = row; } for (unsigned int x = 0; x < IMAGE_WIDTH; x++) { uint8_t pixel = *(rows[2] +x); // Skip water pixels. if (maskRow[x / (IMAGE_WIDTH / MASK_WIDTH)] > 128) { if (pixel == 255) { int unsigned totalWeight = 0; float neighborTotal = 0; for (int i = -2; i <= 2; i++) { for (int j = -2; j <= 2; j++) { int neighborX = x + i; int neighborY = y + j; if (neighborX < 0 || neighborX >= IMAGE_WIDTH || neighborY < 0 || neighborY >= IMAGE_HEIGHT) { continue; } uint8_t neighbor = *(rows[2 + j] + neighborX); if (neighbor == 255) { continue; } int unsigned weight = 5 - (abs(i) + abs(j)); neighborTotal += weight * neighbor; totalWeight += weight; } } if (totalWeight > 11) { pixel = (uint8_t) (neighborTotal / totalWeight); fixedPixels++; } else { unfixablePixels++; } } } else { waterPixels++; } outputBuffer[x] = pixel; } fwrite(outputBuffer, 1, IMAGE_WIDTH, outputFile); if (ferror(outputFile)) { printf("An error occurred writing row %d on pass %d.\n", y, pass); exit(1); } } fclose(outputFile); fclose(inputFile); fclose(maskFile); printf("Finished pass %d in %f seconds. Fixed pixels: %ld\t Unfixable pixels: %ld\t Water pixels: %ld.\n", pass, (clock() - startTime) / (double) CLOCKS_PER_SEC, fixedPixels, unfixablePixels, waterPixels); } int main(int argc, char* argv[]) { if (argc != 2) { printf("Proper usage: %s \n", argv[0]); } int unsigned passes = strtol(argv[1], NULL, 10); if (passes <= 0) { printf("passes argument must be a positive integer.\n"); exit(1); } findPixelRange(&minPixelValue, &maxPixelValue, &bitReductionFactor); reduceBitDepth(); for (int unsigned i = 1; i <= passes; i++) { fillMissingPixels(i); } return 0; }