waterfall.c

   1 #include <fcntl.h>
   2 #include <fftw3.h>
   3 #include <getopt.h>
   4 #include <glib.h>
   5 #include <math.h>
   6 #include <gdk/gdkkeysyms.h>
   7 #include <gtk/gtk.h>
   8 #include <stdio.h>
   9 #include <stdlib.h>
  10 #include <string.h>
  11 #include <unistd.h>
  12
  13 #define SAMPLE_RATE 44100
  14 #define N_SLICES 64
  15 #define SAMPLE_SIZE 2
  16
  17 struct waterfall_context {
  18         GtkDrawingArea *drawingarea;
  19         GtkRuler *frequency_ruler;
  20         GtkRuler *time_ruler;
  21         GdkPixbuf *original; /* A 1-to-1 map from spectral data to pixels. */
  22         GdkPixbuf *resized; /* A backing store for the scaled spectrogram. */
  23         struct {
  24                 GdkPixbuf *part1; /* Slices from zero to the slice cursor. */
  25                 GdkPixbuf *part2; /* Slices from the slice cursor to the end. */
  26                 int part1_size;
  27                 int part2_size;
  28         } subpixbufs;
  29         struct {
  30                 /* The size of the spectrogram drawing area. */
  31                 int width;
  32                 int height;
  33         } size;
  34         struct {
  35                 int f_low;
  36                 int f_high;
  37                 int noise_floor;
  38                 double sensitivity;
  39                 double hertz_scale;
  40         } zoom;
  41         char *buf;
  42         size_t buf_size;
  43         size_t buf_index;
  44         size_t integration_samples;
  45         int slice;
  46         int n_slices;
  47         struct {
  48                 double *samples;
  49                 fftw_complex *spectrum;
  50                 fftw_plan p;
  51         } fft;
  52 };
  53
  54 void invalidate_spectrogram(struct waterfall_context *waterfall);
  55
  56 void waterfall_open_output(GtkFileChooser *chooser, gpointer user_data)
  57 {
  58         char *filename = gtk_file_chooser_get_filename(chooser);
  59         printf("Open %s\n", filename);
  60         g_free(filename);
  61
  62         gtk_widget_hide(GTK_WIDGET(chooser));
  63 }
  64
  65 void waterfall_expose_spectrogram(GtkWidget *widget,
  66                                   GdkEventExpose *event,
  67                                   gpointer user_data)
  68 {
  69         GtkDrawingArea *spectrogram_drawingarea = GTK_DRAWING_AREA(widget);
  70         struct waterfall_context *waterfall = user_data;
  71         cairo_t *cr;
  72
  73         if (!waterfall->resized) {
  74                 double part1_size, part2_size;
  75                 int resolution = waterfall->zoom.f_high - waterfall->zoom.f_low + 1;
  76                 double scale_x = (double) waterfall->size.width / waterfall->n_slices;
  77                 double scale_y = (double) waterfall->size.height / resolution;
  78                 double offset_y = -waterfall->zoom.f_low * scale_y;
  79                 int row, rowstride, channels;
  80                 guchar *pixels;
  81
  82                 part1_size = (double) waterfall->subpixbufs.part1_size / waterfall->n_slices * waterfall->size.width;
  83                 part2_size = waterfall->size.width - (int) part1_size;
  84
  85                 /* Create a new backing store, then render into it. */
  86                 waterfall->resized = gdk_pixbuf_new(GDK_COLORSPACE_RGB, FALSE, 8, waterfall->size.width, waterfall->size.height);
  87                 if (waterfall->subpixbufs.part2) {
  88                         /* Render the oldest slices first. */
  89                         gdk_pixbuf_scale(waterfall->subpixbufs.part2, waterfall->resized,
  90                                          0, 0, (int) part2_size, waterfall->size.height,
  91                                          0, offset_y, scale_x, scale_y,
  92                                          GDK_INTERP_NEAREST);
  93                 }
  94                 if (waterfall->subpixbufs.part1) {
  95                         /* Then render the more recent slices. */
  96                         gdk_pixbuf_scale(waterfall->subpixbufs.part1, waterfall->resized,
  97                                          (int) part2_size, 0, (int) part1_size, waterfall->size.height,
  98                                          part2_size, offset_y, scale_x, scale_y,
  99                                          GDK_INTERP_NEAREST);
 100                 }
 101                 pixels = gdk_pixbuf_get_pixels(waterfall->resized);
 102                 rowstride = gdk_pixbuf_get_rowstride(waterfall->resized);
 103                 channels = gdk_pixbuf_get_n_channels(waterfall->resized);
 104                 /* Draw a narrow marker to represent the border between part1 and part2. */
 105                 for (row = 0; row < waterfall->size.height; row++) {
 106                         pixels[row*rowstride + ((int) part2_size % waterfall->size.width) * channels + 1] = 192;
 107                 }
 108         }
 109
 110         cr = gdk_cairo_create(spectrogram_drawingarea->widget.window);
 111
 112         gdk_cairo_set_source_pixbuf(cr, waterfall->resized, 0, 0);
 113         cairo_paint(cr);
 114
 115         cairo_destroy(cr);
 116 }
 117
 118 void waterfall_configure_spectrogram(GtkWidget *widget,
 119                                      GdkEventConfigure *event,
 120                                      gpointer user_data)
 121 {
 122         struct waterfall_context *waterfall = user_data;
 123
 124         printf("configure %d %d,%d+%d,%d send_event=%d on window %p\n",
 125                event->type, event->x, event->y, event->width, event->height, event->send_event, event->window);
 126
 127         waterfall->size.width = event->width;
 128         waterfall->size.height = event->height;
 129
 130         invalidate_spectrogram(waterfall);
 131 }
 132
 133 gboolean waterfall_key_press(GtkWidget *widget,
 134                              GdkEventKey *event,
 135                              gpointer user_data)
 136 {
 137         struct waterfall_context *waterfall = user_data;
 138         int resolution = waterfall->zoom.f_high - waterfall->zoom.f_low + 1;
 139         int headroom = waterfall->integration_samples/2 - waterfall->zoom.f_high;
 140         gboolean zoom_changed = FALSE;
 141
 142         printf("key press\n");
 143         switch (event->keyval) {
 144         case GDK_KEY_Up:
 145                 printf("up\n");
 146                 waterfall->zoom.f_low -= MIN(waterfall->zoom.f_low, resolution / 2);
 147                 waterfall->zoom.f_high -= MIN(waterfall->zoom.f_low, resolution / 2);
 148                 zoom_changed = TRUE;
 149                 break;
 150         case GDK_KEY_Down:
 151                 printf("down\n");
 152                 waterfall->zoom.f_low += MIN(headroom, resolution / 2);
 153                 waterfall->zoom.f_high += MIN(headroom, resolution / 2);
 154                 zoom_changed = TRUE;
 155                 break;
 156         case GDK_KEY_plus:
 157                 printf("plus\n");
 158                 if (resolution / 4 > 1) {
 159                         waterfall->zoom.f_low += resolution / 4;
 160                         waterfall->zoom.f_high -= resolution / 4;
 161                         zoom_changed = TRUE;
 162                 }
 163                 break;
 164         case GDK_KEY_minus:
 165                 printf("minus\n");
 166                 waterfall->zoom.f_low -= MIN(waterfall->zoom.f_low, resolution / 4);
 167                 waterfall->zoom.f_high += MIN(headroom, resolution / 4);
 168                 zoom_changed = TRUE;
 169                 break;
 170         case GDK_KEY_bracketleft:
 171                 printf("[\n");
 172                 waterfall->zoom.noise_floor--;
 173                 break;
 174         case GDK_KEY_bracketright:
 175                 printf("]\n");
 176                 waterfall->zoom.noise_floor++;
 177                 break;
 178         case GDK_KEY_braceleft:
 179                 printf("{\n");
 180                 waterfall->zoom.sensitivity /= 1.5;
 181                 break;
 182         case GDK_KEY_braceright:
 183                 printf("}\n");
 184                 waterfall->zoom.sensitivity *= 1.5;
 185                 break;
 186         }
 187
 188         if (zoom_changed) {
 189                 invalidate_spectrogram(waterfall);
 190         }
 191
 192         return FALSE;
 193 }
 194
 195 void connect_signal(GtkBuilder *builder,
 196                     GObject *object,
 197                     const gchar *signal_name,
 198                     const gchar *handler_name,
 199                     GObject *connect_object,
 200                     GConnectFlags flags,
 201                     gpointer user_data)
 202 {
 203         static struct {
 204                 char const *name;
 205                 GCallback fn;
 206         } handlers[] = {
 207                 { "gtk_main_quit", G_CALLBACK(&gtk_main_quit) },
 208                 { "gtk_widget_hide_on_delete", G_CALLBACK(&gtk_widget_hide_on_delete) },
 209                 { "gtk_widget_show", G_CALLBACK(&gtk_widget_show) },
 210                 { "waterfall_open_output", G_CALLBACK(&waterfall_open_output) },
 211                 { NULL, NULL }
 212         };
 213         int i;
 214
 215         for (i = 0; handlers[i].name; i++) {
 216                 if (!strcmp(handlers[i].name, handler_name)) {
 217                         GCallback handler = handlers[i].fn;
 218                         g_signal_connect_object(object, signal_name, handler,
 219                                                 connect_object, flags);
 220                         return;
 221                 }
 222         }
 223
 224         /* TODO: Connect to some error-spewing function? */
 225 }
 226
 227 void invalidate_spectrogram(struct waterfall_context *waterfall)
 228 {
 229         GdkRegion *visible_region;
 230
 231         /* Just invalidate everything. */
 232         if (waterfall->resized) {
 233                 g_object_unref(G_OBJECT(waterfall->resized));
 234                 waterfall->resized = NULL;
 235         }
 236         visible_region = gdk_drawable_get_visible_region(waterfall->drawingarea->widget.window);
 237         gdk_window_invalidate_region(waterfall->drawingarea->widget.window,
 238                                      visible_region,
 239                                      FALSE);
 240         gdk_region_destroy(visible_region);
 241
 242         gtk_ruler_set_range(waterfall->frequency_ruler,
 243                             waterfall->zoom.f_low * waterfall->zoom.hertz_scale,
 244                             waterfall->zoom.f_high * waterfall->zoom.hertz_scale,
 245                             waterfall->zoom.f_high * waterfall->zoom.hertz_scale,
 246                             (waterfall->zoom.f_high - waterfall->zoom.f_low) * waterfall->zoom.hertz_scale);
 247 }
 248
 249 void split_pixbuf(struct waterfall_context *waterfall)
 250 {
 251         int max_resolution = waterfall->integration_samples/2 + 1;
 252
 253         if (waterfall->subpixbufs.part1) {
 254                 g_object_unref(G_OBJECT(waterfall->subpixbufs.part1));
 255         }
 256         if (waterfall->subpixbufs.part2) {
 257                 g_object_unref(G_OBJECT(waterfall->subpixbufs.part2));
 258         }
 259
 260         waterfall->subpixbufs.part1_size = waterfall->slice % waterfall->n_slices;
 261         waterfall->subpixbufs.part2_size = waterfall->n_slices - waterfall->subpixbufs.part1_size;
 262         if (waterfall->subpixbufs.part1_size) {
 263                 waterfall->subpixbufs.part1 = gdk_pixbuf_new_subpixbuf(waterfall->original,
 264                                                                        0, 0,
 265                                                                        waterfall->subpixbufs.part1_size, max_resolution);
 266         } else {
 267                 waterfall->subpixbufs.part1 = NULL;
 268         }
 269         if (waterfall->subpixbufs.part2_size) {
 270                 waterfall->subpixbufs.part2 = gdk_pixbuf_new_subpixbuf(waterfall->original,
 271                                                                        waterfall->subpixbufs.part1_size, 0,
 272                                                                        waterfall->subpixbufs.part2_size, max_resolution);
 273         } else {
 274                 waterfall->subpixbufs.part2 = NULL;
 275         }
 276 }
 277
 278 gboolean waterfall_input(GIOChannel *source,
 279                          GIOCondition condition,
 280                          gpointer userdata)
 281 {
 282         struct waterfall_context *waterfall = userdata;
 283         double normalization = log10(waterfall->integration_samples * 32768);
 284         gsize n;
 285         int batch_size = 0;
 286         GError *err = NULL;
 287
 288         switch (g_io_channel_read_chars(source,
 289                                         waterfall->buf,
 290                                         waterfall->buf_size - waterfall->buf_index,
 291                                         &n,
 292                                         &err)) {
 293         case G_IO_STATUS_NORMAL:
 294                 waterfall->buf_index += n;
 295                 break;
 296         case G_IO_STATUS_AGAIN:
 297                 printf("Try again later\n");
 298                 g_error_free(err);
 299                 return TRUE;
 300         case G_IO_STATUS_EOF:
 301                 printf("End of input\n");
 302                 exit(0);
 303                 break;
 304         case G_IO_STATUS_ERROR:
 305                 printf("Error - %s\n", err->message);
 306                 exit(1);
 307                 break;
 308         default:
 309                 exit(1);
 310                 break;
 311         }
 312
 313         while (waterfall->buf_index >= waterfall->integration_samples*SAMPLE_SIZE) {
 314                 int i, row, rowstride, channels, modslice;
 315                 int resolution = waterfall->integration_samples/2 + 1;
 316                 guchar *pixels;
 317
 318                 /* Analyze the spectrum. */
 319                 printf("Integrate %d %zu %d (%d)\n", batch_size, waterfall->buf_index,
 320                        waterfall->slice, waterfall->slice % waterfall->n_slices);
 321                 for (i = 0; i < waterfall->integration_samples; i++) {
 322                         int16_t x;
 323                         memcpy(&x, waterfall->buf + i*SAMPLE_SIZE, sizeof (x));
 324                         waterfall->fft.samples[i] = x;
 325                 }
 326                 fftw_execute(waterfall->fft.p);
 327
 328                 /* Draw a new slice. */
 329                 modslice = waterfall->slice % waterfall->n_slices;
 330                 pixels = gdk_pixbuf_get_pixels(waterfall->original);
 331                 rowstride = gdk_pixbuf_get_rowstride(waterfall->original);
 332                 channels = gdk_pixbuf_get_n_channels(waterfall->original);
 333                 for (row = 0; row < resolution; row++) {
 334                         double signal_i = waterfall->fft.spectrum[row][0], signal_q = waterfall->fft.spectrum[row][1];
 335                         /* XXX Give invsqrt a chance to happen. */
 336                         double power = -log10(1.0 / hypot(signal_i, signal_q)) - normalization;
 337                         double z = (power + waterfall->zoom.noise_floor) * waterfall->zoom.sensitivity;
 338                         if (row > 50 && row < 60) {
 339                                 printf("Power = %g (%g)\n", power, z);
 340                         }
 341                         /* Drawing to the original-size pixbuf, resampling happens in the expose event. */
 342                         pixels[row*rowstride + modslice * channels + 0] = z >= 256 ? MIN(z - 256, 255) : 0;
 343                         pixels[row*rowstride + modslice * channels + 1] = z >= 256 ? MAX(511 - z, 0) : MAX(z, 0);
 344                         pixels[row*rowstride + modslice * channels + 2] = z < 256 ? MIN(255 - z, 255) : 0;
 345                 }
 346                 waterfall->slice++;
 347
 348                 /* Consume input. */
 349                 memmove(waterfall->buf,
 350                         waterfall->buf + waterfall->integration_samples*SAMPLE_SIZE,
 351                         waterfall->buf_size - waterfall->integration_samples*SAMPLE_SIZE);
 352                 waterfall->buf_index -= waterfall->integration_samples*SAMPLE_SIZE;
 353                 batch_size++;
 354         }
 355
 356         /* Scroll the waterfall. */
 357         split_pixbuf(waterfall);
 358         invalidate_spectrogram(waterfall);
 359
 360         if (batch_size > 2) {
 361                 printf("I can't keep up (%d)\n", batch_size);
 362         }
 363
 364         if (err) {
 365                 g_error_free(err);
 366         }
 367
 368         return TRUE;
 369 }
 370
 371 int main(int argc, char *argv[])
 372 {
 373         GtkBuilder *gtk_builder;
 374         GtkWindow *waterfall_window;
 375         GIOChannel *stdin_channel;
 376         int stdin_flags;
 377         GError *err = NULL;
 378         struct waterfall_context waterfall = {
 379                 .buf_index = 0,
 380                 .slice = 0,
 381                 .zoom = {
 382                         .f_low = 0,
 383                         .noise_floor = 5,
 384                         .sensitivity = 200,
 385                 },
 386         };
 387
 388         waterfall.integration_samples = SAMPLE_RATE/10;
 389         waterfall.n_slices = N_SLICES;
 390
 391         while (1) {
 392                 int c = getopt(argc, argv, "H:t:");
 393                 if (c == -1) {
 394                         break;
 395                 }
 396
 397                 switch (c) {
 398                 case 'H':
 399                         waterfall.n_slices = atoi(optarg);
 400                         break;
 401                 case 't':
 402                         waterfall.integration_samples = atof(optarg) * SAMPLE_RATE;
 403                         break;
 404                 default:
 405                         fprintf(stderr, "%s: '%c': Invalid option\n", argv[0], c);
 406                         break;
 407                 }
 408         }
 409
 410         waterfall.zoom.f_high = waterfall.integration_samples/2;
 411         waterfall.zoom.hertz_scale = 21.050 / waterfall.zoom.f_high;
 412
 413         gtk_init(&argc, &argv);
 414
 415         gtk_builder = gtk_builder_new();
 416         gtk_builder_add_from_file(gtk_builder, "waterfall.glade", NULL);
 417         gtk_builder_connect_signals_full(gtk_builder, &connect_signal, NULL);
 418
 419         /* Give signal handlers a means of accessing gtk_builder. */
 420         g_object_set(gtk_builder_get_object(gtk_builder,
 421                                             "waterfall_window"),
 422                      "user-data", gtk_builder,
 423                      NULL);
 424
 425         waterfall_window = GTK_WINDOW(gtk_builder_get_object(gtk_builder, "waterfall_window"));
 426         waterfall.drawingarea = GTK_DRAWING_AREA(gtk_builder_get_object(gtk_builder, "image_detail_drawingarea"));
 427         waterfall.frequency_ruler = GTK_RULER(gtk_builder_get_object(gtk_builder, "frequency_ruler"));
 428         waterfall.time_ruler = GTK_RULER(gtk_builder_get_object(gtk_builder, "time_ruler"));
 429
 430         gtk_ruler_set_range(waterfall.time_ruler,
 431                             N_SLICES * waterfall.integration_samples / 44100.0,
 432                             0,
 433                             N_SLICES * waterfall.integration_samples / 44100.0,
 434                             N_SLICES * waterfall.integration_samples / 44100.0);
 435
 436         g_signal_connect(G_OBJECT(waterfall.drawingarea), "expose-event",
 437                          G_CALLBACK(&waterfall_expose_spectrogram), &waterfall);
 438         g_signal_connect(G_OBJECT(waterfall.drawingarea), "configure-event",
 439                          G_CALLBACK(&waterfall_configure_spectrogram), &waterfall);
 440         g_signal_connect(G_OBJECT(waterfall_window), "key-press-event",
 441                          G_CALLBACK(&waterfall_key_press), &waterfall);
 442
 443         /* Don't block on stdin. */
 444         stdin_flags = fcntl(0, F_GETFL, 0);
 445         fcntl(0, F_SETFL, stdin_flags | O_NONBLOCK);
 446
 447         stdin_channel = g_io_channel_unix_new(0);
 448         g_io_channel_set_encoding(stdin_channel, NULL, &err);
 449         g_io_channel_set_buffer_size(stdin_channel, waterfall.integration_samples * SAMPLE_SIZE);
 450         waterfall.buf = malloc(waterfall.integration_samples * SAMPLE_SIZE);
 451         waterfall.buf_size = waterfall.integration_samples * SAMPLE_SIZE;
 452         g_io_add_watch(stdin_channel, G_IO_IN, &waterfall_input, &waterfall);
 453
 454         waterfall.original = gdk_pixbuf_new(GDK_COLORSPACE_RGB, FALSE, 8,
 455                                             waterfall.n_slices, waterfall.integration_samples/2 + 1);
 456         split_pixbuf(&waterfall);
 457
 458         waterfall.fft.samples = fftw_malloc(sizeof (*waterfall.fft.samples) * waterfall.integration_samples);
 459         waterfall.fft.spectrum = fftw_malloc(sizeof (*waterfall.fft.spectrum) * (waterfall.integration_samples/2 + 1));
 460         waterfall.fft.p = fftw_plan_dft_r2c_1d(waterfall.integration_samples,
 461                                                waterfall.fft.samples, waterfall.fft.spectrum,
 462                                                FFTW_ESTIMATE);
 463
 464         gtk_widget_show(GTK_WIDGET(waterfall_window));
 465
 466         invalidate_spectrogram(&waterfall);
 467
 468         gtk_main();
 469
 470         fftw_destroy_plan(waterfall.fft.p);
 471         fftw_free(waterfall.fft.spectrum);
 472         fftw_free(waterfall.fft.samples);
 473 }