find_coincidences_5byte.c

#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <string.h>
// #include <sys/stat.h>
#include <time.h>

#define TIMESTAMP_BYTE_COUNT 4
#define INITIALIZATION_OFFSET 40

// off_t fsize(const char*);
// int write_results_to_file(const char*, uint8_t(*)[4], uint16_t);

size_t fsize(FILE* file)
{
	if( !file ) { return 0; }

	const size_t curr_position = ftell(file);

	fseek(file, 0, SEEK_END);
	size_t file_size = ftell(file);
	fseek(file, curr_position, SEEK_SET);

	return file_size;
}

int write_results_to_file( const char* file_name,\
						   uint8_t (*results)[4], \
						   uint16_t number_of_elements
						  )
{
	FILE *file;
	file = fopen(file_name, "a");

	for(int i = 0; i < number_of_elements; i++)
	{
		printf("fourfold coincidence: ");
		for(int j = 0; j < 3; j++)
			printf("%d,", results[i][j]/32 + 1);
		printf("%d\n", results[i][3]/32 + 1);

		fprintf(file, "%d,%d,%d,%d\n",\
			results[i][0]/32 + 1,\
			results[i][1]/32 + 1,\
			results[i][2]/32 + 1,\
			results[i][3]/32 + 1
		);
	}

	fclose(file);
	return 0;
}

int main(int argc, char* argv[]) {

	clock_t begin_time;
	begin_time = clock();

	const char DEFAULT_INPUT_FILENAME[] = "timetags5byteLong.bin";
	const char DEFAULT_OUTPUT_FILENAME[] = "processed_timestamps.txt";
	char INPUT_FILENAME[128];
	char OUTPUT_FILENAME[128];
	// if a single input argument is given...
	if( argc == 2 )
		strcpy(INPUT_FILENAME, argv[1]);
	else
		strcpy(INPUT_FILENAME, DEFAULT_INPUT_FILENAME);
	strcpy(OUTPUT_FILENAME, DEFAULT_OUTPUT_FILENAME);

	const uint8_t WINDOW_SIZE = 100;
	// we use arrays with 4 elements, being interested only in fourfold coincidences. 
	uint32_t timestamps[TIMESTAMP_BYTE_COUNT];
	uint8_t channels[TIMESTAMP_BYTE_COUNT];
	// assuming that there are no more than 1000 coincidences
	uint8_t quadruples[1000][TIMESTAMP_BYTE_COUNT];
	uint16_t quadruples_pos = 0;
	// curr_byte and curr_timestamp will be used in the for loop
	uint32_t curr_byte;
	uint32_t curr_timestamp;
	
	// open file in binary read mode
	FILE *file;
	file = fopen(INPUT_FILENAME, "rb");
	if( file == NULL )
	{
		perror(INPUT_FILENAME);
		return(-1);
	}
	// Initialize array of chars able to contain all the bytes.
	// The first 40 bytes are to be ignored.
	// If the number of bytes is not a multiple of 5, the exceeding bytes are
	// not read.
	// In this way, file_size effectively equals 5 times the number of 5-byte blocks
	size_t file_size = fsize(file) - INITIALIZATION_OFFSET;
	file_size = file_size - file_size % 5;
	uint8_t* data = malloc(sizeof(uint8_t) * file_size);
	if (!data) {
		fprintf(stderr,
			"Couldn't allocate %u bytes of data.\n", sizeof(uint8_t) * file_size
		  );
		return 1;
	}

	// ignore the first 40 bytes
	fseek(file, INITIALIZATION_OFFSET, SEEK_SET);
	// the rest of the file is stored in data
	fread(data, sizeof(uint8_t), file_size, file);

	// We start by reading the first 5 bytes stored in data.
	// This is achieved reinterpreting the first 4 bytes in data as a single
	// 4-byte unsigned integer number.
	// The fifth byte encodes the channel.
	timestamps[0] = *(uint32_t*) &data[0];
	channels[0] = data[TIMESTAMP_BYTE_COUNT];
	uint8_t timestamps_pos = 1;

	// Loop through data, at blocks of 5 bytes (remembering that we already
	// processed the first 5 bytes)
	for (int timetag_index = 1;
		 timetag_index < file_size/5;
		 timetag_index++
		)
	{
		// curr_byte points to the first byte of the 5-byte block that is
		// being processed
		curr_byte = timetag_index * 5;
		curr_timestamp = *(uint32_t*) &data[curr_byte];
		if(	curr_timestamp - timestamps[0] <= WINDOW_SIZE )
		{
			timestamps[timestamps_pos] = curr_timestamp;
			channels[timestamps_pos] = data[curr_byte + TIMESTAMP_BYTE_COUNT];
			timestamps_pos++;
		}
		else // we reinitialize the arrays and start a new window
		{
			// first the quadruple is stored, if we have one
			if(timestamps_pos == 4)
			{
				memcpy(&quadruples[quadruples_pos], channels, sizeof(channels));
				// printf("Timetag index = %u\n", timetag_index);
				quadruples_pos++;
			}
			// then a new window is opened
			timestamps_pos = 1;
			timestamps[0] = curr_timestamp;
			channels[0] = data[curr_byte + TIMESTAMP_BYTE_COUNT];
		}
	}
	
	fclose(file);
	free(data);

	printf("Time required: %f\n", (double) (clock()-begin_time) / CLOCKS_PER_SEC);

	write_results_to_file(OUTPUT_FILENAME, quadruples, quadruples_pos);
	
	return 0;

}