some random stuff

.gitignore

@@ -4,5 +4,6 @@
 
 *.o
 *.d
+*.png
 
 /lost

Makefile

@@ -26,9 +26,14 @@ OBJS := $(patsubst %.c,%.o,$(SRCS))
 DEPS := $(patsubst %.c,%.d,$(SRCS))
 BIN  := lost
 
+BSD  := bright-star-database.tsv
+
 LDFLAGS := -lcairo -lm
 
-all: $(BIN)
+all: $(BIN) $(BSD)
+
+$(BSD): download-bright-star-database.sh
+	./download-bright-star-database.sh
 
 $(BIN): $(OBJS)
 	$(CC) $(LDFLAGS) -o $(BIN) $(OBJS)
@@ -40,5 +45,6 @@ $(BIN): $(OBJS)
 
 clean:
 	rm -f $(OBJS) $(DEPS)
+	rm -i $(BSD)
 
 .PHONY: all clean

compile_flags.txt

@@ -0,0 +1 @@
+-std=c99

download-astrometry.sh

@@ -0,0 +1,30 @@
+#!/bin/sh
+
+# Downloads an image and metadata from astrometry, to check astrometry's (supposedly good) results
+# against our (supposedly crappy) ones.
+
+set -e
+
+if [ -z $1 ]
+then
+    echo 'USAGE: ./download-astrometry.sh <user-image-id>'
+    exit 1
+fi
+
+if ! command -v convert >/dev/null
+then
+    echo 'Please install imagemagick for jpeg->png conversion!'
+    exit 1
+fi
+
+astrometry_host=http://nova.astrometry.net
+
+user_image_page=$(curl ""$astrometry_host/user_images/$1)
+axy_link=$astrometry_host$(echo "$user_image_page" | grep -oP "/axy_file/[0-9]+")
+jpeg_link=$astrometry_host$(echo "$user_image_page" | grep -oP '(?<="original": ")/image/[0-9]+')
+
+mkdir "$1"
+curl -o "$1/axy.fits" "$axy_link"
+curl -o "$1/img.jpg" "$jpeg_link"
+convert "$1/img.jpg" "$1/img.png"
+

src/centroiders.c

@@ -2,26 +2,28 @@
 
 #include <stdio.h>
 #include <stdlib.h>
+#include <math.h>
 
 // DUMMY
 
 typedef struct dummy_centroid_config_t {
     int num_stars;
 } dummy_centroid_config_t;
 
-void *pv_dummy_centroid(unsigned char *pc_image,
-                        int i_image_width,
-                        int i_image_height,
-                        int *pi_result_length,
-                        void *pv_config) {
+star_t *pv_dummy_centroid(unsigned char *pc_image,
+                              int i_image_width,
+                              int i_image_height,
+                              int *pi_result_length,
+                              void *pv_config) {
     dummy_centroid_config_t *px_config = (dummy_centroid_config_t *)pv_config;
-    centroid_base_t *px_result = malloc(sizeof(centroid_base_t) * px_config->num_stars);
+    star_t *px_result = calloc(sizeof(star_t), px_config->num_stars);
 
     *pi_result_length = px_config->num_stars;
 
     for (int i = 0; i < px_config->num_stars; i++) {
-        px_result[i].l_x = rand() % i_image_width * 1000000;
-        px_result[i].l_y = rand() % i_image_height * 1000000;
+        px_result[i].f_x = rand() % i_image_width;
+        px_result[i].f_y = rand() % i_image_height;
+        px_result[i].f_radius_x = 10.0;
     }
 
     return (void *)px_result;
@@ -54,3 +56,84 @@ centroid_algorithm_t x_centroid_algorithm(int i) {
         return x_empty_centroid_algorithm;
     }
 }
+
+// OTHER CENTROID RELATED FUNCTIONS
+
+float f_centroid_distance(star_t x_one, star_t x_two) {
+    float f_dist_x = x_one.f_x - x_two.f_x;
+    float f_dist_y = x_one.f_y - x_two.f_y;
+    return sqrt(f_dist_x*f_dist_x + f_dist_y*f_dist_y);
+}
+
+// helper for x_compare_centroids
+static void v_closest_centroids(float f_threshold, star_t *px_one,
+                                int i_one_length, star_t *px_two,
+                                int i_two_length, int *result) {
+    for (int i = 0; i < i_one_length; i++) {
+        float f_distance_closest = INFINITY;
+        int   i_index_closest = -1;
+
+        for (int k = 0; k < i_two_length; k++) {
+            float f_distance_curr = f_centroid_distance(px_one[i], px_two[k]);
+            if (f_distance_curr < f_threshold && f_distance_curr < f_distance_closest) {
+                f_distance_closest = f_distance_curr;
+                i_index_closest = k;
+            }
+        }
+
+        result[i] = i_index_closest;
+    }
+}
+
+centroid_comparison_t x_compare_centroids(float f_threshold,
+                                          star_t *px_expected,
+                                          int i_expected_length,
+                                          star_t *px_actual,
+                                          int i_actual_length) {
+
+    centroid_comparison_t result = { 0 };
+    // maps from indexes in each list to the closest centroid from other list
+    int *pi_expected_to_actual = malloc(sizeof(int) * i_expected_length),
+        *pi_actual_to_expected = malloc(sizeof(int) * i_actual_length);
+
+    v_closest_centroids(f_threshold,
+                        px_expected, i_expected_length,
+                        px_actual, i_actual_length,
+                        pi_expected_to_actual);
+
+    v_closest_centroids(f_threshold,
+                        px_actual, i_actual_length,
+                        px_expected, i_expected_length,
+                        pi_actual_to_expected);
+
+    // any expected stars whose closest actual star does not refer to them are missing
+    for (int i = 0; i < i_expected_length; i++) {
+        if (pi_expected_to_actual[i] == -1 ||
+            pi_actual_to_expected[pi_expected_to_actual[i]] != i) {
+            result.i_stars_missing_num++;
+        } else {
+            result.f_mean_error += f_centroid_distance(px_expected[i],
+                                                       px_actual[pi_expected_to_actual[i]]);
+        }
+    }
+    result.f_mean_error /= (i_expected_length - result.i_stars_missing_num);
+
+
+    // any actual star whose closest expected star does not refer to them is extra
+    for (int i = 0; i < i_actual_length; i++) {
+        if (pi_actual_to_expected[i] == -1 ||
+            pi_expected_to_actual[pi_actual_to_expected[i]] != i) {
+            result.i_stars_extra_num++;
+        }
+    }
+
+    free(pi_expected_to_actual);
+    free(pi_actual_to_expected);
+}
+
+void v_print_centroid_comparison(centroid_comparison_t x_comparison) {
+    printf("Extra stars: %d\nMissing stars: %d\nMean error: %f\n",
+           x_comparison.i_stars_extra_num,
+           x_comparison.i_stars_missing_num,
+           x_comparison.f_mean_error);
+}

src/centroiders.h

@@ -1,37 +1,50 @@
 #ifndef CENTROID_H
 #define CENTROID_H
 
-// centroiding algorithms might return additional information, which only certain algorithms can take advantage of. Those centroiding algorithms should return 
-typedef struct centroid_base_t {
-    long l_x; // pixels*10^6
-    long l_y; // pixels*10^6
-} centroid_base_t;
-
-typedef struct centroid_magnitude_t {
-    centroid_base_t x_base;
-    long l_magnitude;
-} centroid_magnitude_t;
-
-enum {
-    CENTROID_TYPE_BASE,
-    CENTROID_TYPE_MAGNITUDE,
-};
-
-typedef void *(*centroid_function_t)(unsigned char *pc_image,
-                                     int  i_image_width,
-                                     int  i_image_height,
-                                     int  *pi_result_length,
-                                     void *pv_config);
+// Different centroiding algorithms may return different amounts of data. Because there won't ever
+// be /that/ many different centroids in a single image, it's acceptable to waste some space on
+// un-used fields. Set unused fields to negative
+typedef struct star_t {
+    float f_x; // pixels*10^6
+    float f_y; // pixels*10^6
+    long  l_magnitude; // some relative number
+    float f_radius_x;
+    float f_radius_y;  // if omitted, but x is present, assume circular.
+    // eccentricity?
+} star_t;
+
+typedef struct centroid_comparison_t {
+    int i_stars_extra_num;    // stars in actual but not expected. Ideally 0
+    int i_stars_missing_num;  // stars is expected but not actual. Ideally 0
+    float f_mean_error;       // average distance from actual to expected star
+    // I would add 99th percentile or something, but the really far away stars should really just
+    // count in extra_num
+} centroid_comparison_t;
+
+typedef star_t *(*centroid_function_t)(unsigned char *pc_image,
+                                           int  i_image_width,
+                                           int  i_image_height,
+                                           int  *pi_result_length,
+                                           void *pv_config);
 typedef void *(*centroid_config_function_t)(void);
 
 typedef struct centroid_algorithm_t {
     char *s_name;
-    int  i_centroid_type;
 
     centroid_function_t        pf_algorithm;
     centroid_config_function_t pf_config;
 } centroid_algorithm_t;
 
 centroid_algorithm_t x_centroid_algorithm(int i);
 
+float f_centroid_distance(star_t x_one, star_t x_two);
+// compare two lists of centroids to find differences.
+centroid_comparison_t x_compare_centroids(float f_distance_threshold, // stars further apart than
+                                                                      // are different stars.
+                                          star_t *px_expected,
+                                          int i_expected_length,
+                                          star_t *px_actual,
+                                          int i_actual_length);
+void v_print_centroid_comparison(centroid_comparison_t x_comparison);
+
 #endif

src/io.c

@@ -0,0 +1,116 @@
+#include "io.h"
+
+#include <stdio.h>
+#include <inttypes.h>
+#include <math.h>
+
+catalog_t *px_bsd_parse(char *s_path) {
+    catalog_t      *px_result;
+    catalog_star_t *px_star_current;
+    FILE           *px_file;
+    long           l_raj2000_h, l_raj2000_l, // high and low parts
+                   l_dej2000_h, l_dej2000_l;
+    int            i_magnitude_h, i_magnitude_l;
+    char           c_weird;
+
+    px_result = calloc(sizeof(catalog_t), 1);
+
+    px_file = fopen(s_path, "r");
+    if (px_file == NULL) {
+        printf("Error opening file: %s\n", strerror(errno));
+        return NULL;
+    }
+
+    while (EOF != fscanf(px_file, "%ld.%ld|%ld.%ld|%*d|%c|%d.%d",
+                         &l_raj2000_h, &l_raj2000_l,
+                         &l_dej2000_h, &l_dej2000_l,
+                         &c_weird,
+                         &i_magnitude_h, &i_magnitude_l)) {
+
+        if (px_result->l_stars_length % 1000 == 0) {
+            px_result->px_stars = realloc(px_result->px_stars,
+                                          sizeof(catalog_star_t) *
+                                          (px_result->l_stars_length/1000 + 1) * 1000);
+        }
+
+        px_star_current = &px_result->px_stars[px_result->l_stars_length];
+        px_star_current->l_raj2000 = l_raj2000_h * 1000000 + l_raj2000_l;
+        px_star_current->l_dej2000 = l_dej2000_h * 1000000 + l_dej2000_l;
+        px_star_current->i_magnitude = i_magnitude_h * 100 + i_magnitude_l;
+        px_star_current->i_weird = c_weird != ' ';
+        px_result->l_stars_length++;
+    }
+
+    fclose(px_file);
+    return px_result;
+}
+
+unsigned char *pc_surface_grayscale(cairo_surface_t *px_surface) {
+    int i_width, i_height;
+    unsigned char *pc_result;
+    uint32_t      *px_surface_data, x_pixel;
+
+    if (cairo_image_surface_get_format(px_surface) != CAIRO_FORMAT_ARGB32 &&
+        cairo_image_surface_get_format(px_surface) != CAIRO_FORMAT_RGB24) {
+        puts("Can't convert weird image formats to grayscale.");
+        return NULL;
+    }
+
+    i_width  = cairo_image_surface_get_width(px_surface);
+    i_height = cairo_image_surface_get_height(px_surface);
+
+    pc_result = malloc(i_width * i_height);
+    px_surface_data = (uint32_t *)cairo_image_surface_get_data(px_surface);
+
+    for (int i = 0; i < i_height * i_width; i++) {
+        x_pixel = px_surface_data[i];
+        // use "luminosity" method of grayscaling
+        pc_result[i] = round(
+            (x_pixel>>16 &0xFF) *0.21 +
+            (x_pixel>>8  &0xFF) *0.71 +
+            (x_pixel     &0xFF) *0.07
+            );
+    }
+
+    return pc_result;
+}
+
+void v_surface_plot_centroids(cairo_surface_t *px_surface,
+                              star_t *px_centroids,
+                              int i_centroids_length,
+                              double red,
+                              double green,
+                              double blue,
+                              double alpha) {
+    cairo_t *px_cairo;
+
+    px_cairo = cairo_create(px_surface);
+    cairo_set_source_rgba(px_cairo, red, green, blue, alpha);
+    cairo_set_line_width(px_cairo, 1.0); // I wonder what <1.0 does?
+
+    while (i_centroids_length --> 0) {
+        if (px_centroids->f_radius_x > 0.0f) {
+            float f_radius_x = px_centroids->f_radius_x;
+            float f_radius_y = px_centroids->f_radius_y > 0.0f ?
+                px_centroids->f_radius_y : f_radius_x;
+
+            // Rectangles should be entirely /outside/ the radius of the star, so the star is fully
+            // visible.
+            cairo_rectangle(px_cairo,
+                            floor(px_centroids->f_x - f_radius_x) - 1,
+                            floor(px_centroids->f_y - f_radius_y) - 1,
+                            ceil(f_radius_x) * 2 + 2,
+                            ceil(f_radius_y) * 2 + 2);
+            cairo_stroke(px_cairo);
+        } else {
+            cairo_rectangle(px_cairo,
+                            floor(px_centroids->f_x),
+                            floor(px_centroids->f_y),
+                            1, 1);
+            cairo_fill(px_cairo);
+        }
+        px_centroids++;
+    }
+    cairo_destroy(px_cairo);
+}
+