Browse Source

Switch to tabs

tags/v1.8.0
Xerbo 4 years ago
parent
commit
a35a4c602b
10 changed files with 486 additions and 486 deletions
  1. +2
    -2
      Makefile
  2. +85
    -85
      color.c
  3. +75
    -75
      dsp.c
  4. +24
    -24
      filter.c
  5. +66
    -66
      image.c
  6. +39
    -39
      main.c
  7. +15
    -15
      median.c
  8. +44
    -44
      pngio.c
  9. +48
    -48
      reg.c
  10. +88
    -88
      satcal.h

+ 2
- 2
Makefile View File

@@ -9,9 +9,9 @@ aptdec: $(OBJS)

reg.o: reg.c
color.o: color.c
main.o: main.c
main.o: main.c
media.o: median.c
dsp.o: dsp.c
dsp.o: dsp.c
filter.o: filter.c
image.o: image.c
pngio.o: pngio.c


+ 85
- 85
color.c View File

@@ -30,96 +30,96 @@ rgb_t applyPalette(char *palette, int val){
}

rgb_t RGBcomposite(rgb_t top, float top_a, rgb_t bottom, float bottom_a){
return (rgb_t){
MCOMPOSITE(top.r, top_a, bottom.r, bottom_a),
MCOMPOSITE(top.g, top_a, bottom.g, bottom_a),
MCOMPOSITE(top.b, top_a, bottom.b, bottom_a)
};
return (rgb_t){
MCOMPOSITE(top.r, top_a, bottom.r, bottom_a),
MCOMPOSITE(top.g, top_a, bottom.g, bottom_a),
MCOMPOSITE(top.b, top_a, bottom.b, bottom_a)
};
}

// The "I totally didn't just steal this from WXtoImg" palette
char TempPalette[256*3] = {
"\x45\x0\x8f\x46\x0\x91\x47\x0\x92\x48\x0\x94\x49\x0\x96\x4a\x0\x98\x4b\x0\x9b\x4d\x0\x9d"
"\x4e\x0\xa0\x50\x0\xa2\x51\x0\xa5\x52\x0\xa7\x54\x0\xaa\x56\x0\xae\x57\x0\xb1"
"\x58\x0\xb4\x5a\x0\xb7\x5c\x0\xba\x5e\x0\xbd\x5f\x0\xc0\x61\x0\xc4\x64\x0\xc8"
"\x66\x0\xcb\x68\x0\xce\x69\x0\xd1\x68\x0\xd4\x65\x0\xd7\x63\x0\xda\x61\x0\xdd"
"\x5d\x0\xe1\x5b\x0\xe4\x59\x0\xe6\x56\x0\xe9\x53\x0\xeb\x50\x0\xee\x4d\x0\xf0"
"\x49\x0\xf3\x47\x0\xfc\x43\x0\xfa\x31\x0\xbf\x20\x0\x89\x20\x0\x92\x1e\x0\x95"
"\x1b\x0\x97\x19\x0\x9a\x17\x0\x9c\x15\x0\x9e\x12\x0\xa0\xf\x0\xa3\xf\x2\xa5"
"\xe\x6\xa8\xe\xa\xab\xe\xd\xad\xe\x11\xb1\xd\x15\xb4\xd\x18\xb7\xd\x1c\xba"
"\xb\x21\xbd\xa\x25\xc0\xa\x29\xc3\x9\x2d\xc6\x8\x33\xca\x7\x36\xcd\x7\x3b\xd0"
"\x7\x41\xd3\x5\x45\xd6\x4\x4b\xd9\x4\x50\xdc\x3\x55\xde\x2\x5d\xe2\x1\x61\xe5"
"\x0\x66\xe7\x0\x6c\xea\x0\x72\xec\x0\x78\xee\x0\x7d\xf0\x0\x82\xf3\x0\x8d\xfc"
"\x0\x90\xfa\x0\x71\xbf\x0\x54\x89\x0\x5c\x91\x0\x61\x94\x0\x64\x96\x0\x68\x97"
"\x0\x6d\x99\x0\x71\x9b\x0\x75\x9d\x0\x79\x9f\x0\x7e\xa0\x0\x82\xa2\x0\x87\xa4"
"\x0\x8c\xa6\x0\x92\xaa\x0\x96\xac\x0\x9c\xae\x0\xa2\xb1\x0\xa6\xb3\x0\xaa\xb5"
"\x0\xad\xb7\x0\xb1\xba\x0\xb6\xbe\x0\xba\xc0\x0\xbe\xc2\x0\xc2\xc5\x0\xc6\xc6"
"\x0\xca\xc9\x0\xcc\xca\x0\xcf\xcb\x0\xd2\xcc\x0\xd4\xcc\x0\xd6\xcc\x0\xd9\xcb"
"\x0\xdb\xcb\x0\xde\xcb\x0\xe0\xcb\x0\xe2\xcc\x0\xea\xd2\x0\xea\xcf\x0\xb9\xa4"
"\x0\x8e\x7a\x1\x94\x7c\x4\x97\x79\x7\x99\x75\x9\x9b\x71\xd\x9d\x6b\x10\x9f\x67"
"\x12\xa1\x63\x15\xa3\x5f\x17\xa5\x59\x1a\xa8\x55\x1d\xaa\x50\x20\xac\x4b\x24\xaf\x45"
"\x28\xb2\x41\x2b\xb5\x3b\x2e\xb8\x35\x31\xba\x30\x34\xbd\x2b\x39\xbf\x24\x3f\xc1\x17"
"\x49\xc5\x8\x4f\xc8\x1\x4f\xca\x0\x4e\xcd\x0\x4e\xcf\x0\x4f\xd2\x0\x54\xd5\x0"
"\x5d\xd8\x0\x68\xdb\x0\x6e\xdd\x0\x74\xdf\x0\x7a\xe2\x0\x7f\xe4\x0\x85\xe7\x0"
"\x8b\xe9\x0\x8f\xeb\x0\x9b\xf3\x0\x9e\xf2\x0\x7e\xbb\x0\x60\x8a\x0\x68\x92\x0"
"\x6d\x95\x0\x71\x96\x0\x75\x98\x0\x7b\x9a\x0\x7f\x9d\x0\x83\x9f\x0\x87\xa1\x0"
"\x8c\xa2\x0\x8f\xa5\x0\x92\xa7\x0\x96\xa9\x0\x9a\xad\x0\x9d\xb0\x0\xa1\xb2\x0"
"\xa5\xb5\x0\xa9\xb7\x0\xad\xba\x0\xb2\xbd\x0\xb6\xbf\x0\xbb\xc3\x0\xbf\xc6\x0"
"\xc3\xc8\x0\xc8\xcb\x0\xcc\xce\x0\xd0\xd1\x0\xd3\xd2\x0\xd5\xd4\x0\xd9\xd4\x0"
"\xdc\xd4\x0\xde\xd5\x0\xe1\xd5\x0\xe3\xd5\x0\xe6\xd4\x0\xe8\xd1\x0\xea\xce\x0"
"\xf2\xcf\x0\xf2\xca\x0\xbb\x99\x0\x8a\x6e\x0\x92\x72\x0\x95\x72\x0\x97\x71\x0"
"\x9a\x70\x0\x9c\x6e\x0\x9e\x6d\x0\xa0\x6b\x0\xa3\x6a\x0\xa5\x68\x0\xa8\x67\x0"
"\xab\x66\x0\xae\x65\x0\xb2\x63\x0\xb4\x61\x0\xb7\x5f\x0\xba\x5d\x0\xbd\x5c\x0"
"\xc0\x59\x0\xc3\x57\x0\xc6\x54\x0\xca\x50\x0\xcd\x4d\x0\xd0\x4a\x0\xd3\x47\x0"
"\xd6\x43\x0\xd9\x40\x0\xdc\x3d\x0\xde\x39\x0\xe2\x33\x0\xe5\x2f\x0\xe7\x2c\x0"
"\xea\x28\x0\xec\x23\x0\xef\x1f\x0\xf1\x1a\x0\xf3\x14\x0\xfb\xf\x0\xfa\xd\x0"
"\xc1\x5\x0\x8e\x0\x0\x97\x0\x0\x9b\x0\x0\x9e\x0\x0\xa1\x0\x0\xa5\x0\x0"
"\xa9\x0\x0\xad\x0\x0\xb1\x0\x0\xb6\x0\x0\xba\x0\x0\xbd\x0\x0\xc2\x0\x0"
"\xc8\x0\x0\xcc\x0\x0\xcc\x0\x0"
"\x45\x0\x8f\x46\x0\x91\x47\x0\x92\x48\x0\x94\x49\x0\x96\x4a\x0\x98\x4b\x0\x9b\x4d\x0\x9d"
"\x4e\x0\xa0\x50\x0\xa2\x51\x0\xa5\x52\x0\xa7\x54\x0\xaa\x56\x0\xae\x57\x0\xb1"
"\x58\x0\xb4\x5a\x0\xb7\x5c\x0\xba\x5e\x0\xbd\x5f\x0\xc0\x61\x0\xc4\x64\x0\xc8"
"\x66\x0\xcb\x68\x0\xce\x69\x0\xd1\x68\x0\xd4\x65\x0\xd7\x63\x0\xda\x61\x0\xdd"
"\x5d\x0\xe1\x5b\x0\xe4\x59\x0\xe6\x56\x0\xe9\x53\x0\xeb\x50\x0\xee\x4d\x0\xf0"
"\x49\x0\xf3\x47\x0\xfc\x43\x0\xfa\x31\x0\xbf\x20\x0\x89\x20\x0\x92\x1e\x0\x95"
"\x1b\x0\x97\x19\x0\x9a\x17\x0\x9c\x15\x0\x9e\x12\x0\xa0\xf\x0\xa3\xf\x2\xa5"
"\xe\x6\xa8\xe\xa\xab\xe\xd\xad\xe\x11\xb1\xd\x15\xb4\xd\x18\xb7\xd\x1c\xba"
"\xb\x21\xbd\xa\x25\xc0\xa\x29\xc3\x9\x2d\xc6\x8\x33\xca\x7\x36\xcd\x7\x3b\xd0"
"\x7\x41\xd3\x5\x45\xd6\x4\x4b\xd9\x4\x50\xdc\x3\x55\xde\x2\x5d\xe2\x1\x61\xe5"
"\x0\x66\xe7\x0\x6c\xea\x0\x72\xec\x0\x78\xee\x0\x7d\xf0\x0\x82\xf3\x0\x8d\xfc"
"\x0\x90\xfa\x0\x71\xbf\x0\x54\x89\x0\x5c\x91\x0\x61\x94\x0\x64\x96\x0\x68\x97"
"\x0\x6d\x99\x0\x71\x9b\x0\x75\x9d\x0\x79\x9f\x0\x7e\xa0\x0\x82\xa2\x0\x87\xa4"
"\x0\x8c\xa6\x0\x92\xaa\x0\x96\xac\x0\x9c\xae\x0\xa2\xb1\x0\xa6\xb3\x0\xaa\xb5"
"\x0\xad\xb7\x0\xb1\xba\x0\xb6\xbe\x0\xba\xc0\x0\xbe\xc2\x0\xc2\xc5\x0\xc6\xc6"
"\x0\xca\xc9\x0\xcc\xca\x0\xcf\xcb\x0\xd2\xcc\x0\xd4\xcc\x0\xd6\xcc\x0\xd9\xcb"
"\x0\xdb\xcb\x0\xde\xcb\x0\xe0\xcb\x0\xe2\xcc\x0\xea\xd2\x0\xea\xcf\x0\xb9\xa4"
"\x0\x8e\x7a\x1\x94\x7c\x4\x97\x79\x7\x99\x75\x9\x9b\x71\xd\x9d\x6b\x10\x9f\x67"
"\x12\xa1\x63\x15\xa3\x5f\x17\xa5\x59\x1a\xa8\x55\x1d\xaa\x50\x20\xac\x4b\x24\xaf\x45"
"\x28\xb2\x41\x2b\xb5\x3b\x2e\xb8\x35\x31\xba\x30\x34\xbd\x2b\x39\xbf\x24\x3f\xc1\x17"
"\x49\xc5\x8\x4f\xc8\x1\x4f\xca\x0\x4e\xcd\x0\x4e\xcf\x0\x4f\xd2\x0\x54\xd5\x0"
"\x5d\xd8\x0\x68\xdb\x0\x6e\xdd\x0\x74\xdf\x0\x7a\xe2\x0\x7f\xe4\x0\x85\xe7\x0"
"\x8b\xe9\x0\x8f\xeb\x0\x9b\xf3\x0\x9e\xf2\x0\x7e\xbb\x0\x60\x8a\x0\x68\x92\x0"
"\x6d\x95\x0\x71\x96\x0\x75\x98\x0\x7b\x9a\x0\x7f\x9d\x0\x83\x9f\x0\x87\xa1\x0"
"\x8c\xa2\x0\x8f\xa5\x0\x92\xa7\x0\x96\xa9\x0\x9a\xad\x0\x9d\xb0\x0\xa1\xb2\x0"
"\xa5\xb5\x0\xa9\xb7\x0\xad\xba\x0\xb2\xbd\x0\xb6\xbf\x0\xbb\xc3\x0\xbf\xc6\x0"
"\xc3\xc8\x0\xc8\xcb\x0\xcc\xce\x0\xd0\xd1\x0\xd3\xd2\x0\xd5\xd4\x0\xd9\xd4\x0"
"\xdc\xd4\x0\xde\xd5\x0\xe1\xd5\x0\xe3\xd5\x0\xe6\xd4\x0\xe8\xd1\x0\xea\xce\x0"
"\xf2\xcf\x0\xf2\xca\x0\xbb\x99\x0\x8a\x6e\x0\x92\x72\x0\x95\x72\x0\x97\x71\x0"
"\x9a\x70\x0\x9c\x6e\x0\x9e\x6d\x0\xa0\x6b\x0\xa3\x6a\x0\xa5\x68\x0\xa8\x67\x0"
"\xab\x66\x0\xae\x65\x0\xb2\x63\x0\xb4\x61\x0\xb7\x5f\x0\xba\x5d\x0\xbd\x5c\x0"
"\xc0\x59\x0\xc3\x57\x0\xc6\x54\x0\xca\x50\x0\xcd\x4d\x0\xd0\x4a\x0\xd3\x47\x0"
"\xd6\x43\x0\xd9\x40\x0\xdc\x3d\x0\xde\x39\x0\xe2\x33\x0\xe5\x2f\x0\xe7\x2c\x0"
"\xea\x28\x0\xec\x23\x0\xef\x1f\x0\xf1\x1a\x0\xf3\x14\x0\xfb\xf\x0\xfa\xd\x0"
"\xc1\x5\x0\x8e\x0\x0\x97\x0\x0\x9b\x0\x0\x9e\x0\x0\xa1\x0\x0\xa5\x0\x0"
"\xa9\x0\x0\xad\x0\x0\xb1\x0\x0\xb6\x0\x0\xba\x0\x0\xbd\x0\x0\xc2\x0\x0"
"\xc8\x0\x0\xcc\x0\x0\xcc\x0\x0"
};

char PrecipPalette[256*3] = {
"\xe0\x98\x8\xec\x84\x10\xf5\x70\x1b\xfc\x5c\x29\xff\x49\x38\xff\x37\x4a"
"\xfb\x28\x5d\xf5\x1a\x71\xeb\xf\x85\xdf\x8\x99\xd0\x3\xad\xc0\x1\xbf"
"\xad\x3\xd0\x9a\x7\xdf\x86\xf\xeb\x72\x1a\xf5\x5e\x27\xfb\x4b\x37\xfe"
"\x39\x48\xff\x29\x5b\xfc\x1b\x6f\xf5\x10\x83\xec\x8\x98\xe0\x3\xab\xd2"
"\x1\xbe\xc1\x2\xcf\xaf\x7\xdd\x9c\xe\xea\x88\x19\xf4\x73\x26\xfb\x5f"
"\x35\xfe\x4c\x47\xff\x3a\x5a\xfc\x2a\x6d\xf6\x1d\x82\xed\x11\x96\xe1\x9"
"\xa9\xd3\x3\xbc\xc3\x1\xcd\xb1\x2\xdc\x9e\x6\xe9\x89\xe\xf3\x75\x18"
"\xfa\x61\x25\xfe\x4e\x34\xff\x3c\x45\xfc\x2c\x58\xf7\x1e\x6c\xee\x12\x80"
"\xe3\x9\x94\xd5\x4\xa8\xc4\x1\xba\xb2\x2\xcc\x9f\x6\xdb\x8b\xd\xe8"
"\x77\x17\xf2\x63\x23\xfa\x50\x33\xfe\x3d\x44\xff\x2d\x56\xfd\x1f\x6a\xf7"
"\x13\x7e\xef\xa\x92\xe4\x4\xa6\xd6\x1\xb9\xc6\x2\xca\xb4\x5\xda\xa1"
"\xc\xe7\x8d\x16\xf1\x79\x22\xf9\x65\x31\xfe\x51\x42\xff\x3f\x54\xfd\x2f"
"\x68\xf8\x20\x7c\xf0\x14\x90\xe5\xb\xa4\xd7\x4\xb7\xc7\x1\xc9\xb6\x1"
"\xd8\xa3\x5\xe6\x8f\xb\xf1\x7b\x15\xf9\x66\x21\xfd\x53\x30\xff\x41\x40"
"\xfd\x30\x53\xf9\x21\x66\xf1\x15\x7a\xe6\xb\x8f\xd9\x5\xa3\xc9\x1\xb6"
"\xb7\x1\xc7\xa5\x4\xd7\x91\xb\xe5\x7c\x14\xf0\x68\x20\xf8\x55\x2e\xfd"
"\x42\x3f\xff\x31\x51\xfe\x22\x64\xf9\x16\x78\xf2\xc\x8d\xe7\x5\xa1\xda"
"\x2\xb4\xca\x1\xc6\xb9\x4\xd6\xa6\xa\xe4\x93\x13\xef\x7e\x1f\xf7\x6a"
"\x2d\xfd\x56\x3d\xff\x44\x4f\xfe\x33\x63\xfa\x24\x77\xf2\x17\x8b\xe8\xd"
"\x9f\xdb\x6\xb2\xcc\x2\xc4\xbb\x1\xd4\xa8\x4\xe2\x94\x9\xee\x80\x12"
"\xf7\x6c\x1d\xfc\x58\x2b\xff\x45\x3c\xfe\x34\x4e\xfa\x25\x61\xf3\x18\x75"
"\xe9\xe\x89\xdc\x6\x9d\xcd\x2\xb1\xbc\x1\xc3\xaa\x3\xd3\x96\x9\xe1"
"\x82\x11\xed\x6e\x1c\xf6\x5a\x2a\xfc\x47\x3a\xff\x36\x4c\xfe\x26\x5f\xfb"
"\x19\x73\xf4\xf\x87\xea\x7\x9b\xde\x2\xaf\xcf\x1\xc1\xbe\x3\xd2\xac"
"\x8\xe0\x98\x10\xec\x84\x1b\xf5\x6f\x29\xfc\x5c\x39\xff\x49\x4a\xff\x37"
"\x5d\xfb\x28\x71\xf5\x1a\x85\xeb\xf\x9a\xdf\x7\xad\xd0\x3\xbf\xc0\x1"
"\xd0\xad\x3\xdf\x9a\x7\xeb\x86\xf\xf5\x71\x1a\xfb\x5d\x27\xff\x4a\x37"
"\xff\x39\x49\xfc\x29\x5c\xf5\x1b\x6f\xec\x10\x84\xe0\x8\x98\xd2\x3\xab"
"\xc1\x1\xbe\xaf\x2\xcf\x9c\x7\xde\x87\xe\xea\x73\x19\xf4\x5f\x26\xfb"
"\x4c\x36\xfe\x3a\x47\xff\x2a\x5a\xfc\x1c\x6e\xf6\x11\x82\xed\x9\x96\xe1"
"\x3\xaa\xd3\x1\xbc\xc3\x2\xcd\xb1\x6\xdc\x9d\xe\xe9\x89\x18\xf3\x75"
"\x25\xfa\x61\x34\xfe\x4e\x45\xff\x3c\x58\xfc\x2c\x6c\xf7\x1e\x80\xee\x12"
"\x94\xe2\x9\xa8\xd4\x4\xbb\xc4\x1\xcc\xb2\x2\xdb\x9f\x6\xe8\x8b\xd"
"\xf2\x77\x17\xfa\x63\x24\xfe\x4f\x33\xff\x3d\x44\xfd\x2d\x56\xf7\x1f\x6a"
"\xef\x13\x7e\xe4\xa\x92\xd6\x4\xa6\xc6\x1\xb9\xb4\x2\xca\xa1\x5\xda"
"\x8d\xc\xe7\x79\x16\xf2\x64\x22\xf9\x51\x31\xfe\x3f\x42\xff\x2e\x55\xfd"
"\x20\x68\xf8\x14\x7c\xf0\xb\x91\xe5\x4\xa4\xd7\x1\xb7\xc7\x1\xc9\xb6"
"\x5\xd9\xa3\xb\xe6\x8f\x15\xf1\x7a\x21\xf9\x66\x30\xfd\x53\x41\xff\x40"
"\x53\xfd\x30\x66\xf9\x21\x7b\xf1\x15\x8f\xe6\xb\xa3\xd8\x5\xb6\xc9\x1"
"\xc7\xb7\x1\xd7\xa4\x4\xe5\x91\xb\xf0\x7c\x14\xf8\x68\x20\xfd\x54\x2e"
"\xff\x42\x3f\xfe\x31\x51\xf9\x22\x65\xf1\x16\x79\xe7\xc\x8d\xda\x5\xa1"
"\xca\x2\xb4\xb9\x1\xc6\xa6\x4\xd6\x92\xa\xe4\x7e\x13\xef\x6a\x1f\xf7"
"\x56\x2d\xfd\x44\x3d\xff\x33\x4f\xfe\x24\x63\xfa"
"\xe0\x98\x8\xec\x84\x10\xf5\x70\x1b\xfc\x5c\x29\xff\x49\x38\xff\x37\x4a"
"\xfb\x28\x5d\xf5\x1a\x71\xeb\xf\x85\xdf\x8\x99\xd0\x3\xad\xc0\x1\xbf"
"\xad\x3\xd0\x9a\x7\xdf\x86\xf\xeb\x72\x1a\xf5\x5e\x27\xfb\x4b\x37\xfe"
"\x39\x48\xff\x29\x5b\xfc\x1b\x6f\xf5\x10\x83\xec\x8\x98\xe0\x3\xab\xd2"
"\x1\xbe\xc1\x2\xcf\xaf\x7\xdd\x9c\xe\xea\x88\x19\xf4\x73\x26\xfb\x5f"
"\x35\xfe\x4c\x47\xff\x3a\x5a\xfc\x2a\x6d\xf6\x1d\x82\xed\x11\x96\xe1\x9"
"\xa9\xd3\x3\xbc\xc3\x1\xcd\xb1\x2\xdc\x9e\x6\xe9\x89\xe\xf3\x75\x18"
"\xfa\x61\x25\xfe\x4e\x34\xff\x3c\x45\xfc\x2c\x58\xf7\x1e\x6c\xee\x12\x80"
"\xe3\x9\x94\xd5\x4\xa8\xc4\x1\xba\xb2\x2\xcc\x9f\x6\xdb\x8b\xd\xe8"
"\x77\x17\xf2\x63\x23\xfa\x50\x33\xfe\x3d\x44\xff\x2d\x56\xfd\x1f\x6a\xf7"
"\x13\x7e\xef\xa\x92\xe4\x4\xa6\xd6\x1\xb9\xc6\x2\xca\xb4\x5\xda\xa1"
"\xc\xe7\x8d\x16\xf1\x79\x22\xf9\x65\x31\xfe\x51\x42\xff\x3f\x54\xfd\x2f"
"\x68\xf8\x20\x7c\xf0\x14\x90\xe5\xb\xa4\xd7\x4\xb7\xc7\x1\xc9\xb6\x1"
"\xd8\xa3\x5\xe6\x8f\xb\xf1\x7b\x15\xf9\x66\x21\xfd\x53\x30\xff\x41\x40"
"\xfd\x30\x53\xf9\x21\x66\xf1\x15\x7a\xe6\xb\x8f\xd9\x5\xa3\xc9\x1\xb6"
"\xb7\x1\xc7\xa5\x4\xd7\x91\xb\xe5\x7c\x14\xf0\x68\x20\xf8\x55\x2e\xfd"
"\x42\x3f\xff\x31\x51\xfe\x22\x64\xf9\x16\x78\xf2\xc\x8d\xe7\x5\xa1\xda"
"\x2\xb4\xca\x1\xc6\xb9\x4\xd6\xa6\xa\xe4\x93\x13\xef\x7e\x1f\xf7\x6a"
"\x2d\xfd\x56\x3d\xff\x44\x4f\xfe\x33\x63\xfa\x24\x77\xf2\x17\x8b\xe8\xd"
"\x9f\xdb\x6\xb2\xcc\x2\xc4\xbb\x1\xd4\xa8\x4\xe2\x94\x9\xee\x80\x12"
"\xf7\x6c\x1d\xfc\x58\x2b\xff\x45\x3c\xfe\x34\x4e\xfa\x25\x61\xf3\x18\x75"
"\xe9\xe\x89\xdc\x6\x9d\xcd\x2\xb1\xbc\x1\xc3\xaa\x3\xd3\x96\x9\xe1"
"\x82\x11\xed\x6e\x1c\xf6\x5a\x2a\xfc\x47\x3a\xff\x36\x4c\xfe\x26\x5f\xfb"
"\x19\x73\xf4\xf\x87\xea\x7\x9b\xde\x2\xaf\xcf\x1\xc1\xbe\x3\xd2\xac"
"\x8\xe0\x98\x10\xec\x84\x1b\xf5\x6f\x29\xfc\x5c\x39\xff\x49\x4a\xff\x37"
"\x5d\xfb\x28\x71\xf5\x1a\x85\xeb\xf\x9a\xdf\x7\xad\xd0\x3\xbf\xc0\x1"
"\xd0\xad\x3\xdf\x9a\x7\xeb\x86\xf\xf5\x71\x1a\xfb\x5d\x27\xff\x4a\x37"
"\xff\x39\x49\xfc\x29\x5c\xf5\x1b\x6f\xec\x10\x84\xe0\x8\x98\xd2\x3\xab"
"\xc1\x1\xbe\xaf\x2\xcf\x9c\x7\xde\x87\xe\xea\x73\x19\xf4\x5f\x26\xfb"
"\x4c\x36\xfe\x3a\x47\xff\x2a\x5a\xfc\x1c\x6e\xf6\x11\x82\xed\x9\x96\xe1"
"\x3\xaa\xd3\x1\xbc\xc3\x2\xcd\xb1\x6\xdc\x9d\xe\xe9\x89\x18\xf3\x75"
"\x25\xfa\x61\x34\xfe\x4e\x45\xff\x3c\x58\xfc\x2c\x6c\xf7\x1e\x80\xee\x12"
"\x94\xe2\x9\xa8\xd4\x4\xbb\xc4\x1\xcc\xb2\x2\xdb\x9f\x6\xe8\x8b\xd"
"\xf2\x77\x17\xfa\x63\x24\xfe\x4f\x33\xff\x3d\x44\xfd\x2d\x56\xf7\x1f\x6a"
"\xef\x13\x7e\xe4\xa\x92\xd6\x4\xa6\xc6\x1\xb9\xb4\x2\xca\xa1\x5\xda"
"\x8d\xc\xe7\x79\x16\xf2\x64\x22\xf9\x51\x31\xfe\x3f\x42\xff\x2e\x55\xfd"
"\x20\x68\xf8\x14\x7c\xf0\xb\x91\xe5\x4\xa4\xd7\x1\xb7\xc7\x1\xc9\xb6"
"\x5\xd9\xa3\xb\xe6\x8f\x15\xf1\x7a\x21\xf9\x66\x30\xfd\x53\x41\xff\x40"
"\x53\xfd\x30\x66\xf9\x21\x7b\xf1\x15\x8f\xe6\xb\xa3\xd8\x5\xb6\xc9\x1"
"\xc7\xb7\x1\xd7\xa4\x4\xe5\x91\xb\xf0\x7c\x14\xf8\x68\x20\xfd\x54\x2e"
"\xff\x42\x3f\xfe\x31\x51\xf9\x22\x65\xf1\x16\x79\xe7\xc\x8d\xda\x5\xa1"
"\xca\x2\xb4\xb9\x1\xc6\xa6\x4\xd6\x92\xa\xe4\x7e\x13\xef\x6a\x1f\xf7"
"\x56\x2d\xfd\x44\x3d\xff\x33\x4f\xfe\x24\x63\xfa"
};

+ 75
- 75
dsp.c View File

@@ -74,25 +74,25 @@ static inline double Phase(double I, double Q) {

if(I == 0.0 && Q == 0.0) return 0.0;

if (Q < 0) {
if (Q < 0) {
s = -1;
Q = -Q;
} else {
} else {
s = 1;
}
if (I >= 0) {
r = (I - Q) / (I + Q);
angle = 0.25 - 0.25 * r;
} else {
r = (I + Q) / (Q - I);
angle = 0.75 - 0.25 * r;
}
if(s > 0){
return angle;
}else{
return -angle;
}
if (I >= 0) {
r = (I - Q) / (I + Q);
angle = 0.25 - 0.25 * r;
} else {
r = (I + Q) / (Q - I);
angle = 0.75 - 0.25 * r;
}
if(s > 0){
return angle;
}else{
return -angle;
}
}

@@ -102,44 +102,44 @@ static inline double Phase(double I, double Q) {
*/
static double pll(double I, double Q) {
// PLL coefficient
static double PhaseOsc = 0.0;
double Io, Qo;
double Ip, Qp;
double DPhi;
static double PhaseOsc = 0.0;
double Io, Qo;
double Ip, Qp;
double DPhi;

// Quadrature oscillator / reference
Io = cos(PhaseOsc);
Qo = sin(PhaseOsc);
Io = cos(PhaseOsc);
Qo = sin(PhaseOsc);

// Phase detector
Ip = I * Io + Q * Qo;
Qp = Q * Io - I * Qo;
DPhi = Phase(Ip, Qp);
Ip = I * Io + Q * Qo;
Qp = Q * Io - I * Qo;
DPhi = Phase(Ip, Qp);

// Loop filter
PhaseOsc += 2.0 * M_PI * (K1 * DPhi + FreqOsc);
if (PhaseOsc > M_PI)
PhaseOsc += 2.0 * M_PI * (K1 * DPhi + FreqOsc);
if (PhaseOsc > M_PI)
PhaseOsc -= 2.0 * M_PI;
if (PhaseOsc <= -M_PI)
if (PhaseOsc <= -M_PI)
PhaseOsc += 2.0 * M_PI;

FreqOsc += K2 * DPhi;
if (FreqOsc > ((Fc + DFc) / Fe))
FreqOsc += K2 * DPhi;
if (FreqOsc > ((Fc + DFc) / Fe))
FreqOsc = (Fc + DFc) / Fe;
if (FreqOsc < ((Fc - DFc) / Fe))
if (FreqOsc < ((Fc - DFc) / Fe))
FreqOsc = (Fc - DFc) / Fe;

return Ip;
return Ip;
}

// Convert samples into pixels
static int getamp(double *ampbuff, int count) {
static float inbuff[BLKIN];
static int idxin = 0;
static int nin = 0;
static float inbuff[BLKIN];
static int idxin = 0;
static int nin = 0;

for (int n = 0; n < count; n++) {
double I, Q;
for (int n = 0; n < count; n++) {
double I, Q;

// Get some more samples when needed
if (nin < IQFilterLen * 2 + 2) {
@@ -147,11 +147,11 @@ static int getamp(double *ampbuff, int count) {
int res;
memmove(inbuff, &(inbuff[idxin]), nin * sizeof(float));
idxin = 0;
// Read some samples
res = getsample(&(inbuff[nin]), BLKIN - nin);
nin += res;
// Make sure there is enough samples to continue
if (nin < IQFilterLen * 2 + 2)
return n;
@@ -164,25 +164,25 @@ static int getamp(double *ampbuff, int count) {
// Increment current sample
idxin++;
nin--;
}
}

return count;
return count;
}

// Sub-pixel offsetting + FIR compensation
int getpixelv(float *pvbuff, int count) {
// Amplitude buffer
static double ampbuff[BLKAMP];
static int nam = 0;
static int idxam = 0;
static double ampbuff[BLKAMP];
static int nam = 0;
static int idxam = 0;

double mult;
double mult;

// Gaussian resampling factor
mult = (double) Fi / Fe * FreqLine;
int m = RSFilterLen / mult + 1;
mult = (double) Fi / Fe * FreqLine;
int m = RSFilterLen / mult + 1;

for (int n = 0; n < count; n++) {
for (int n = 0; n < count; n++) {
int shift;

if (nam < m) {
@@ -203,53 +203,53 @@ int getpixelv(float *pvbuff, int count) {

idxam += shift;
nam -= shift;
}
return count;
}
return count;
}

// Get an entire row of pixels, aligned with sync markers
// FIXME: skips noisy lines with no findable sync marker
int getpixelrow(float *pixelv, int nrow, int *zenith) {
static float pixels[PixelLine + SyncFilterLen];
static int npv;
static int synced = 0;
static double max = 0.0;
static float pixels[PixelLine + SyncFilterLen];
static int npv;
static int synced = 0;
static double max = 0.0;
static double minDoppler = 100;

double corr, ecorr, lcorr;
int res;
double corr, ecorr, lcorr;
int res;

// Move the row buffer into the the image buffer
if (npv > 0)
if (npv > 0)
memmove(pixelv, pixels, npv * sizeof(float));

// Get the sync line
if (npv < SyncFilterLen + 2) {
if (npv < SyncFilterLen + 2) {
res = getpixelv(&(pixelv[npv]), SyncFilterLen + 2 - npv);
npv += res;
if (npv < SyncFilterLen + 2)
return 0;
}
}

// Calculate the frequency offset
ecorr = fir(pixelv, Sync, SyncFilterLen);
corr = fir(&pixelv[1], Sync, SyncFilterLen - 1);
lcorr = fir(&pixelv[2], Sync, SyncFilterLen - 2);
FreqLine = 1.0+((ecorr-lcorr) / corr / PixelLine / 4.0);
// Calculate the frequency offset
ecorr = fir(pixelv, Sync, SyncFilterLen);
corr = fir(&pixelv[1], Sync, SyncFilterLen - 1);
lcorr = fir(&pixelv[2], Sync, SyncFilterLen - 2);
FreqLine = 1.0+((ecorr-lcorr) / corr / PixelLine / 4.0);

// Find the point in which ecorr and lcorr intercept
if(fabs(lcorr - ecorr) < minDoppler && fabs(lcorr - ecorr) > 2){
minDoppler = fabs(lcorr - ecorr);
*zenith = nrow;
}
// The point in which the pixel offset is recalculated
if (corr < 0.75 * max) {
if (corr < 0.75 * max) {
synced = 0;
FreqLine = 1.0;
}
max = corr;
}
max = corr;

if (synced < 8) {
int mshift;
@@ -285,20 +285,20 @@ int getpixelrow(float *pixelv, int nrow, int *zenith) {
}

// Get the rest of this row
if (npv < PixelLine) {
if (npv < PixelLine) {
res = getpixelv(&(pixelv[npv]), PixelLine - npv);
npv += res;
if (npv < PixelLine)
return 0;
}
}

// Move the sync lines into the output buffer with the calculated offset
if (npv == PixelLine) {
if (npv == PixelLine) {
npv = 0;
} else {
} else {
memmove(pixels, &(pixelv[PixelLine]), (npv - PixelLine) * sizeof(float));
npv -= PixelLine;
}
}

return 1;
return 1;
}

+ 24
- 24
filter.c View File

@@ -23,45 +23,45 @@

// Finite impulse response
float fir(float *buff, const float *coeff, const int len) {
double r;
double r;

r = 0.0;
for (int i = 0; i < len; i++) {
r += buff[i] * coeff[i];
}
return r;
r = 0.0;
for (int i = 0; i < len; i++) {
r += buff[i] * coeff[i];
}
return r;
}

/* IQ finite impulse response
* Turn samples into a single IQ sample
*/
void iqfir(float *buff, const float *coeff, const int len, double *I, double *Q) {
double i = 0.0, q = 0.0;
double i = 0.0, q = 0.0;

for (int k = 0; k < len; k++) {
q += buff[2*k] * coeff[k];
i += buff[2*k];
}
for (int k = 0; k < len; k++) {
q += buff[2*k] * coeff[k];
i += buff[2*k];
}

i = buff[len-1] - (i / len);
*I = i, *Q = q;
i = buff[len-1] - (i / len);
*I = i, *Q = q;
}

/* Gaussian finite impulse responce compensation
* https://www.recordingblogs.com/wiki/gaussian-window
*/
float rsfir(double *buff, const float *coeff, const int len, const double offset, const double delta) {
double out;
double out;

out = 0.0;
double n = offset;
for (int i = 0; i < (len-1)/delta-1; n += delta, i++) {
int k;
double alpha;
out = 0.0;
double n = offset;
for (int i = 0; i < (len-1)/delta-1; n += delta, i++) {
int k;
double alpha;

k = (int)floor(n);
alpha = n - k;
out += buff[i] * (coeff[k] * (1.0 - alpha) + coeff[k + 1] * alpha);
}
return out;
k = (int)floor(n);
alpha = n - k;
out += buff[i] * (coeff[k] * (1.0 - alpha) + coeff[k + 1] * alpha);
}
return out;
}

+ 66
- 66
image.c View File

@@ -35,22 +35,22 @@ extern void polyreg(const int m, const int n, const double x[], const double y[]

// Compute regression
static void rgcomp(double x[16], rgparam_t * rgpr) {
// { 0.106, 0.215, 0.324, 0.433, 0.542, 0.652, 0.78, 0.87, 0.0 }
const double y[9] = { 31.07, 63.02, 94.96, 126.9, 158.86, 191.1, 228.62, 255.0, 0.0 };
// { 0.106, 0.215, 0.324, 0.433, 0.542, 0.652, 0.78, 0.87, 0.0 }
const double y[9] = { 31.07, 63.02, 94.96, 126.9, 158.86, 191.1, 228.62, 255.0, 0.0 };

polyreg(REGORDER, 9, x, y, rgpr->cf);
polyreg(REGORDER, 9, x, y, rgpr->cf);
}

// Convert a value to 0-255 based off the provided regression curve
static double rgcal(float x, rgparam_t *rgpr) {
double y, p;
double y, p;
int i;
for (i = 0, y = 0.0, p = 1.0; i < REGORDER + 1; i++) {
for (i = 0, y = 0.0, p = 1.0; i < REGORDER + 1; i++) {
y += rgpr->cf[i] * p;
p = p * x;
}
return y;
}
return y;
}

static double tele[16];
@@ -119,38 +119,38 @@ double teleNoise(double wedges[16]){
double noise = 0;
for(int i = 0; i < 9; i++)
noise += fabs(wedges[i] - pattern[i]);
return noise;
}

// Get telemetry data for thermal calibration/equalization
int calibrate(float **prow, int nrow, int offset, int width) {
double teleline[MAX_HEIGHT] = { 0.0 };
double wedge[16];
rgparam_t regr[MAX_HEIGHT/FRAME_LEN + 1];
int telestart, mtelestart = 0;
double wedge[16];
rgparam_t regr[MAX_HEIGHT/FRAME_LEN + 1];
int telestart, mtelestart = 0;
int channel = -1;

// The minimum rows required to decode a full frame
if (nrow < 192) {
if (nrow < 192) {
fprintf(stderr, "Telemetry decoding error, not enough rows\n");
return 0;
}
}

// Calculate average of a row of telemetry
for (int y = 0; y < nrow; y++) {
for (int y = 0; y < nrow; y++) {
for (int x = 3; x < 43; x++)
teleline[y] += prow[y][x + offset + width];
teleline[y] /= 40.0;
}
}

/* Wedge 7 is white and 8 is black, this will have the largest
* difference in brightness, this can be used to find the current
* position within the telemetry.
*/
float max = 0.0f;
for (int n = nrow / 3 - 64; n < 2 * nrow / 3 - 64; n++) {
for (int n = nrow / 3 - 64; n < 2 * nrow / 3 - 64; n++) {
float df;

// (sum 4px below) / (sum 4px above)
@@ -162,20 +162,20 @@ int calibrate(float **prow, int nrow, int offset, int width) {
mtelestart = n;
max = df;
}
}
}

telestart = (mtelestart - 64) % FRAME_LEN;
telestart = (mtelestart - 64) % FRAME_LEN;

// Make sure that theres at least one full frame in the image
if (nrow < telestart + FRAME_LEN) {
// Make sure that theres at least one full frame in the image
if (nrow < telestart + FRAME_LEN) {
fprintf(stderr, "Telemetry decoding error, not enough rows\n");
return 0;
}
}

// Find the least noisy frame
double minNoise = -1;
int bestFrame = -1;
for (int n = telestart, k = 0; n < nrow - FRAME_LEN; n += FRAME_LEN, k++) {
for (int n = telestart, k = 0; n < nrow - FRAME_LEN; n += FRAME_LEN, k++) {
int j;

for (j = 0; j < 16; j++) {
@@ -228,7 +228,7 @@ int calibrate(float **prow, int nrow, int offset, int width) {
}
Cs = rgcal(Cs / i, &regr[k]);
}
}
}

if(bestFrame == -1){
fprintf(stderr, "Something has gone very wrong, please file a bug report.");
@@ -237,8 +237,8 @@ int calibrate(float **prow, int nrow, int offset, int width) {

calibrateImage(prow, nrow, offset, width, regr[bestFrame]);

return channel + 1;
return channel + 1;
}

void distrib(options_t *opts, image_t *img, char *chid) {
@@ -261,7 +261,7 @@ void distrib(options_t *opts, image_t *img, char *chid) {

for(int n = 0; n < img->nrow; n++) {
float *pixelv = img->prow[n];
for(int i = 0; i < CH_WIDTH; i++) {
int y = CLIP((int)pixelv[i + CHA_OFFSET], 0, 255);
int x = CLIP((int)pixelv[i + CHB_OFFSET], 0, 255);
@@ -275,7 +275,7 @@ void distrib(options_t *opts, image_t *img, char *chid) {
for(int x = 0; x < 256; x++)
for(int y = 0; y < 256; y++)
distrib.prow[y][x] = distrib.prow[y][x] / max * 255.0;
extern int ImageOut(options_t *opts, image_t *img, int offset, int width, char *desc, char *chid, char *palette);
ImageOut(&options, &distrib, 0, 256, "Distribution", chid, NULL);
}
@@ -319,21 +319,21 @@ void flipImage(image_t *img, int width, int offset){
#include "satcal.h"

typedef struct {
double Nbb;
double Cs;
double Cb;
int ch;
double Nbb;
double Cs;
double Cb;
int ch;
} tempparam_t;

// IR channel temperature compensation
static void tempcomp(double t[16], int ch, int satnum, tempparam_t *tpr) {
double Tbb, T[4];
double C;
double Tbb, T[4];
double C;

tpr->ch = ch - 4;
tpr->ch = ch - 4;

// Compute equivalent T blackbody temperature
for (int n = 0; n < 4; n++) {
// Compute equivalent T blackbody temperature
for (int n = 0; n < 4; n++) {
float d0, d1, d2;

C = t[9 + n] * 4.0;
@@ -344,57 +344,57 @@ static void tempcomp(double t[16], int ch, int satnum, tempparam_t *tpr) {
T[n] += d1 * C;
C *= C;
T[n] += d2 * C;
}
Tbb = (T[0] + T[1] + T[2] + T[3]) / 4.0;
Tbb = satcal[satnum].rad[tpr->ch].A + satcal[satnum].rad[tpr->ch].B * Tbb;
}
Tbb = (T[0] + T[1] + T[2] + T[3]) / 4.0;
Tbb = satcal[satnum].rad[tpr->ch].A + satcal[satnum].rad[tpr->ch].B * Tbb;

// Compute blackbody radiance temperature
C = satcal[satnum].rad[tpr->ch].vc;
tpr->Nbb = c1 * C * C * C / (expm1(c2 * C / Tbb));
// Compute blackbody radiance temperature
C = satcal[satnum].rad[tpr->ch].vc;
tpr->Nbb = c1 * C * C * C / (expm1(c2 * C / Tbb));

// Store blackbody count and space
tpr->Cs = Cs * 4.0;
tpr->Cb = t[14] * 4.0;
// Store blackbody count and space
tpr->Cs = Cs * 4.0;
tpr->Cb = t[14] * 4.0;
}

// IR channel temperature calibration
static double tempcal(float Ce, int satnum, tempparam_t * rgpr) {
double Nl, Nc, Ns, Ne;
double T, vc;
double Nl, Nc, Ns, Ne;
double T, vc;

Ns = satcal[satnum].cor[rgpr->ch].Ns;
Nl = Ns + (rgpr->Nbb - Ns) * (rgpr->Cs - Ce * 4.0) / (rgpr->Cs - rgpr->Cb);
Nc = satcal[satnum].cor[rgpr->ch].b[0] +
Ns = satcal[satnum].cor[rgpr->ch].Ns;
Nl = Ns + (rgpr->Nbb - Ns) * (rgpr->Cs - Ce * 4.0) / (rgpr->Cs - rgpr->Cb);
Nc = satcal[satnum].cor[rgpr->ch].b[0] +
satcal[satnum].cor[rgpr->ch].b[1] * Nl +
satcal[satnum].cor[rgpr->ch].b[2] * Nl * Nl;

Ne = Nl + Nc;
Ne = Nl + Nc;

vc = satcal[satnum].rad[rgpr->ch].vc;
T = c2 * vc / log(c1 * vc * vc * vc / Ne + 1.0);
T = (T - satcal[satnum].rad[rgpr->ch].A) / satcal[satnum].rad[rgpr->ch].B;
vc = satcal[satnum].rad[rgpr->ch].vc;
T = c2 * vc / log(c1 * vc * vc * vc / Ne + 1.0);
T = (T - satcal[satnum].rad[rgpr->ch].A) / satcal[satnum].rad[rgpr->ch].B;

// Convert to celsius
T -= 273.15;
// Rescale to 0-255 for -100°C to +60°C
T = (T + 100.0) / 160.0 * 255.0;
// Rescale to 0-255 for -100°C to +60°C
T = (T + 100.0) / 160.0 * 255.0;

return T;
}

// Temperature calibration wrapper
void temperature(options_t *opts, image_t *img, int offset, int width){
tempparam_t temp;
tempparam_t temp;

printf("Temperature... ");
fflush(stdout);
printf("Temperature... ");
fflush(stdout);

tempcomp(tele, img->chB, opts->satnum - 15, &temp);
tempcomp(tele, img->chB, opts->satnum - 15, &temp);

for (int y = 0; y < img->nrow; y++) {
for (int x = 0; x < width; x++) {
for (int y = 0; y < img->nrow; y++) {
for (int x = 0; x < width; x++) {
img->prow[y][x + offset] = tempcal(img->prow[y][x + offset], opts->satnum - 15, &temp);
}
}
printf("Done\n");
}
printf("Done\n");
}

+ 39
- 39
main.c View File

@@ -66,7 +66,7 @@ static int processAudio(char *filename, options_t *opts);
static void usage(void);

int main(int argc, char **argv) {
fprintf(stderr, VERSION"\n");
fprintf(stderr, VERSION"\n");

// Check if there are actually any input files
if(argc == optind || argc == 1){
@@ -78,7 +78,7 @@ int main(int argc, char **argv) {

// Parse arguments
int opt;
while ((opt = getopt(argc, argv, "m:d:i:s:e:r")) != EOF) {
while ((opt = getopt(argc, argv, "m:d:i:s:e:r")) != EOF) {
switch (opt) {
case 'd':
opts.path = optarg;
@@ -105,14 +105,14 @@ int main(int argc, char **argv) {
default:
usage();
}
}
}

// Process the files
for (; optind < argc; optind++) {
processAudio(argv[optind], &opts);
}

exit(0);
exit(0);
}

static int processAudio(char *filename, options_t *opts){
@@ -124,7 +124,7 @@ static int processAudio(char *filename, options_t *opts){
char *id[7];
char *name[7];
} ch = {
{ "?", "1", "2", "3A", "4", "5", "3B" },
{ "?", "1", "2", "3A", "4", "5", "3B" },
{ "unknown", "visble", "near-infrared", "mid-infrared", "thermal-infrared", "thermal-infrared", "mid-infrared" }
};

@@ -163,7 +163,7 @@ static int processAudio(char *filename, options_t *opts){
for (img.nrow = 0; img.nrow < MAX_HEIGHT; img.nrow++) {
// Allocate memory for this row
img.prow[img.nrow] = (float *) malloc(sizeof(float) * 2150);
// Write into memory and break the loop when there are no more samples to read
if (getpixelrow(img.prow[img.nrow], img.nrow, &zenith) == 0)
break;
@@ -259,40 +259,40 @@ static int processAudio(char *filename, options_t *opts){
// Distribution image
if (CONTAINS(opts->type, 'd'))
distrib(opts, &img, "d");
return 1;
}

static int initsnd(char *filename) {
SF_INFO infwav;
int res;
SF_INFO infwav;
int res;

// Open audio file
infwav.format = 0;
audioFile = sf_open(filename, SFM_READ, &infwav);
if (audioFile == NULL) {
infwav.format = 0;
audioFile = sf_open(filename, SFM_READ, &infwav);
if (audioFile == NULL) {
fprintf(stderr, "Could not open %s for reading\n", filename);
return 0;
}
}

res = init_dsp(infwav.samplerate);
res = init_dsp(infwav.samplerate);
printf("Input file: %s\n", filename);
if(res < 0) {
if(res < 0) {
fprintf(stderr, "Input sample rate too low: %d\n", infwav.samplerate);
return 0;
}else if(res > 0) {
}else if(res > 0) {
fprintf(stderr, "Input sample rate too high: %d\n", infwav.samplerate);
return 0;
}
printf("Input sample rate: %d\n", infwav.samplerate);
}
printf("Input sample rate: %d\n", infwav.samplerate);

// TODO: accept stereo audio
if (infwav.channels != 1) {
if (infwav.channels != 1) {
fprintf(stderr, "Too many channels in input file: %d\n", infwav.channels);
return 0;
}
}

return 1;
return 1;
}

// Read samples from the wave file
@@ -301,28 +301,28 @@ int getsample(float *sample, int nb) {
}

static void usage(void) {
fprintf(stderr,
fprintf(stderr,
"Aptdec [options] audio files ...\n"
"Options:\n"
" -e [t|h|d|p|f|l] Effects\n"
" t: Crop telemetry\n"
" h: Histogram equalise\n"
" d: Denoise\n"
" p: Precipitation\n"
" f: Flip image\n"
" l: Linear equalise\n"
" t: Crop telemetry\n"
" h: Histogram equalise\n"
" d: Denoise\n"
" p: Precipitation\n"
" f: Flip image\n"
" l: Linear equalise\n"
" -i [r|a|b|c|t|m] Output image\n"
" r: Raw\n"
" a: Channel A\n"
" b: Channel B\n"
" c: False color\n"
" t: Temperature\n"
" m: MCIR\n"
" -d <dir> Image destination directory.\n"
" -s [15-19] Satellite number\n"
" -m <file> Map file\n"
" -r Realtime decode\n"
" r: Raw\n"
" a: Channel A\n"
" b: Channel B\n"
" c: False color\n"
" t: Temperature\n"
" m: MCIR\n"
" -d <dir> Image destination directory.\n"
" -s [15-19] Satellite number\n"
" -m <file> Map file\n"
" -r Realtime decode\n"
"\nRefer to the README for more infomation\n");

exit(EINVAL);
exit(EINVAL);
}

+ 15
- 15
median.c View File

@@ -7,24 +7,24 @@

#define ELEM_SWAP(a, b) { register float t = (a); (a) = (b); (b) = t; }
float quick_select(float arr[], int n) {
int low, median, high;
int middle, ll, hh;
low = 0; high = n-1; median = (low + high) / 2;
for (;;) {
if (high <= low) /* One element only */
return arr[median] ;
if (high == low + 1) { /* Two elements only */
if (arr[low] > arr[high])
ELEM_SWAP(arr[low], arr[high]);
return arr[median];
}
int low, median, high;
int middle, ll, hh;
low = 0; high = n-1; median = (low + high) / 2;
for (;;) {
if (high <= low) /* One element only */
return arr[median] ;
if (high == low + 1) { /* Two elements only */
if (arr[low] > arr[high])
ELEM_SWAP(arr[low], arr[high]);
return arr[median];
}

/* Find median of low, middle and high items; swap into position low */
middle = (low + high) / 2;
if (arr[middle] > arr[high]) ELEM_SWAP(arr[middle], arr[high]);
if (arr[low] > arr[high]) ELEM_SWAP(arr[low], arr[high]);
if (arr[low] > arr[high]) ELEM_SWAP(arr[low], arr[high]);
if (arr[middle] > arr[low]) ELEM_SWAP(arr[middle], arr[low]);

/* Swap low item (now in position middle) into position (low+1) */
@@ -51,6 +51,6 @@ float quick_select(float arr[], int n) {
low = ll;
if (hh >= median)
high = hh - 1;
}
}
}
#undef ELEM_SWAP

+ 44
- 44
pngio.c View File

@@ -191,8 +191,8 @@ png_text meta[] = {

int ImageOut(options_t *opts, image_t *img, int offset, int width, char *desc, char *chid, char *palette){
char outName[384];
sprintf(outName, "%s/%s-%s.png", opts->path, img->name, chid);
sprintf(outName, "%s/%s-%s.png", opts->path, img->name, chid);
meta[1].text = desc;
meta[1].text_length = sizeof(desc);

@@ -208,24 +208,24 @@ int ImageOut(options_t *opts, image_t *img, int offset, int width, char *desc, c
}

// Create writer
png_structp png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);
if (!png_ptr) {
png_structp png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);
if (!png_ptr) {
png_destroy_write_struct(&png_ptr, (png_infopp) NULL);
fprintf(stderr, "Could not create a PNG writer\n");
return 0;
}
png_infop info_ptr = png_create_info_struct(png_ptr);
if (!info_ptr) {
}
png_infop info_ptr = png_create_info_struct(png_ptr);
if (!info_ptr) {
png_destroy_write_struct(&png_ptr, (png_infopp) NULL);
fprintf(stderr, "Could not create a PNG writer\n");
return 0;
}
}

if(palette == NULL && !CONTAINS(opts->effects, 'p') && !fc && opts->map[0] == '\0' && strcmp(chid, "MCIR") != 0){
greyscale = 1;

// Greyscale image
png_set_IHDR(png_ptr, info_ptr, width, img->nrow,
png_set_IHDR(png_ptr, info_ptr, width, img->nrow,
8, PNG_COLOR_TYPE_GRAY, PNG_INTERLACE_NONE,
PNG_COMPRESSION_TYPE_DEFAULT, PNG_FILTER_TYPE_DEFAULT);
}else{
@@ -235,17 +235,17 @@ int ImageOut(options_t *opts, image_t *img, int offset, int width, char *desc, c
PNG_COMPRESSION_TYPE_DEFAULT, PNG_FILTER_TYPE_DEFAULT);
}

png_set_text(png_ptr, info_ptr, meta, 3);
png_set_pHYs(png_ptr, info_ptr, 3636, 3636, PNG_RESOLUTION_METER);
png_set_text(png_ptr, info_ptr, meta, 3);
png_set_pHYs(png_ptr, info_ptr, 3636, 3636, PNG_RESOLUTION_METER);

// Init I/O
pngfile = fopen(outName, "wb");
if (!pngfile) {
pngfile = fopen(outName, "wb");
if (!pngfile) {
fprintf(stderr, "Could not open %s for writing\n", outName);
return 1;
}
png_init_io(png_ptr, pngfile);
png_write_info(png_ptr, info_ptr);
}
png_init_io(png_ptr, pngfile);
png_write_info(png_ptr, info_ptr);

// Move prow into crow, crow ~ color rows
rgb_t *crow[img->nrow];
@@ -287,10 +287,10 @@ int ImageOut(options_t *opts, image_t *img, int offset, int width, char *desc, c
printf("Writing %s", outName);

// Build image
for (int y = 0; y < img->nrow; y++) {
for (int y = 0; y < img->nrow; y++) {
png_color pix[width]; // Color
png_byte mpix[width]; // Mono
int skip = 0;
for (int x = 0; x < width; x++) {
if(skiptele){
@@ -320,21 +320,21 @@ int ImageOut(options_t *opts, image_t *img, int offset, int width, char *desc, c
};
}
}
if(greyscale){
png_write_row(png_ptr, (png_bytep) mpix);
}else{
png_write_row(png_ptr, (png_bytep) pix);
}
}
}

// Tidy up
png_write_end(png_ptr, info_ptr);
fclose(pngfile);
printf("\nDone\n");
png_destroy_write_struct(&png_ptr, &info_ptr);
png_write_end(png_ptr, info_ptr);
fclose(pngfile);
printf("\nDone\n");
png_destroy_write_struct(&png_ptr, &info_ptr);

return 1;
return 1;
}

// TODO: clean up everthing below this comment
@@ -344,44 +344,44 @@ FILE *rt_pngfile;

int initWriter(options_t *opts, image_t *img, int width, int height, char *desc, char *chid){
char outName[384];
sprintf(outName, "%s/%s-%s.png", opts->path, img->name, chid);
sprintf(outName, "%s/%s-%s.png", opts->path, img->name, chid);

meta[1].text = desc;
meta[1].text_length = sizeof(desc);

// Create writer
rt_png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);
if (!rt_png_ptr) {
rt_png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);
if (!rt_png_ptr) {
png_destroy_write_struct(&rt_png_ptr, (png_infopp) NULL);
fprintf(stderr, "Could not create a PNG writer\n");
return 0;
}
rt_info_ptr = png_create_info_struct(rt_png_ptr);
if (!rt_info_ptr) {
}
rt_info_ptr = png_create_info_struct(rt_png_ptr);
if (!rt_info_ptr) {
png_destroy_write_struct(&rt_png_ptr, (png_infopp) NULL);
fprintf(stderr, "Could not create a PNG writer\n");
return 0;
}
}

// Greyscale image
png_set_IHDR(rt_png_ptr, rt_info_ptr, width, height,
png_set_IHDR(rt_png_ptr, rt_info_ptr, width, height,
8, PNG_COLOR_TYPE_GRAY, PNG_INTERLACE_NONE,
PNG_COMPRESSION_TYPE_DEFAULT, PNG_FILTER_TYPE_DEFAULT);
png_set_text(rt_png_ptr, rt_info_ptr, meta, 3);

// Channel = 25cm wide
png_set_pHYs(rt_png_ptr, rt_info_ptr, 3636, 3636, PNG_RESOLUTION_METER);
png_set_pHYs(rt_png_ptr, rt_info_ptr, 3636, 3636, PNG_RESOLUTION_METER);

// Init I/O
rt_pngfile = fopen(outName, "wb");
if (!rt_pngfile) {
rt_pngfile = fopen(outName, "wb");
if (!rt_pngfile) {
fprintf(stderr, "Could not open %s for writing\n", outName);
return 0;
}
png_init_io(rt_png_ptr, rt_pngfile);
png_write_info(rt_png_ptr, rt_info_ptr);
}
png_init_io(rt_png_ptr, rt_pngfile);
png_write_info(rt_png_ptr, rt_info_ptr);
return 1;
}

@@ -389,12 +389,12 @@ void pushRow(float *row, int width){
png_byte pix[width];
for(int i = 0; i < width; i++)
pix[i] = row[i];
png_write_row(rt_png_ptr, (png_bytep) pix);
}

void closeWriter(){
png_write_end(rt_png_ptr, rt_info_ptr);
fclose(rt_pngfile);
png_destroy_write_struct(&rt_png_ptr, &rt_info_ptr);
fclose(rt_pngfile);
png_destroy_write_struct(&rt_png_ptr, &rt_info_ptr);
}

+ 48
- 48
reg.c View File

@@ -18,18 +18,18 @@
static void FactPiv(int N, double A[DMAX][DMAX], double B[], double Cf[]);

void polyreg(const int M, const int N, const double X[], const double Y[], double C[]) {
int R, K, J; /* Loop counters */
double A[DMAX][DMAX]; /* A */
double B[DMAX];
double P[2 * DMAX + 1];
double x, y;
double p;
/* Zero the array */
for (R = 0; R < M + 1; R++)
int R, K, J; /* Loop counters */
double A[DMAX][DMAX]; /* A */
double B[DMAX];
double P[2 * DMAX + 1];
double x, y;
double p;
/* Zero the array */
for (R = 0; R < M + 1; R++)
B[R] = 0;

/* Compute the column vector */
/* Compute the column vector */
for (K = 0; K < N; K++) {
y = Y[K];
x = X[K];
@@ -38,50 +38,50 @@ void polyreg(const int M, const int N, const double X[], const double Y[], doubl
B[R] += y * p;
p = p * x;
}
}
}

/* Zero the array */
for (J = 1; J <= 2 * M; J++)
/* Zero the array */
for (J = 1; J <= 2 * M; J++)
P[J] = 0;
P[0] = N;

/* Compute the sum of powers of x_(K-1) */
for (K = 0; K < N; K++) {
P[0] = N;

/* Compute the sum of powers of x_(K-1) */
for (K = 0; K < N; K++) {
x = X[K];
p = X[K];
for (J = 1; J <= 2 * M; J++) {
P[J] += p;
p = p * x;
}
}
}

/* Determine the matrix entries */
for (R = 0; R < M + 1; R++) {
/* Determine the matrix entries */
for (R = 0; R < M + 1; R++) {
for (K = 0; K < M + 1; K++)
A[R][K] = P[R + K];
}
}

/* Solve the linear system of M + 1 equations: A*C = B
for the coefficient vector C = (c_1,c_2,..,c_M,c_(M+1)) */
FactPiv(M + 1, A, B, C);
/* Solve the linear system of M + 1 equations: A*C = B
for the coefficient vector C = (c_1,c_2,..,c_M,c_(M+1)) */
FactPiv(M + 1, A, B, C);
} /* end main */


/*--------------------------------------------------------*/
static void FactPiv(int N, double A[DMAX][DMAX], double B[], double Cf[]) {
int K, P, C, J; /* Loop counters */
int Row[NMAX]; /* Field with row-number */
double X[DMAX], Y[DMAX];
double SUM, DET = 1.0;
int T;
/* Initialize the pointer vector */
for (J = 0; J < NMAX; J++)
int K, P, C, J; /* Loop counters */
int Row[NMAX]; /* Field with row-number */
double X[DMAX], Y[DMAX];
double SUM, DET = 1.0;
int T;
/* Initialize the pointer vector */
for (J = 0; J < NMAX; J++)
Row[J] = J;

/* Start LU factorization */
for (P = 0; P < N - 1; P++) {
/* Start LU factorization */
for (P = 0; P < N - 1; P++) {
/* Find pivot element */
for (K = P + 1; K < N; K++) {
if (fabs(A[Row[K]][P]) > fabs(A[Row[P]][P])) {
@@ -109,34 +109,34 @@ static void FactPiv(int N, double A[DMAX][DMAX], double B[], double Cf[]) {
A[Row[K]][C] -= A[Row[K]][P] * A[Row[P]][C];
}
}
}
}
/* End of L*U factorization routine */
DET = DET * A[Row[N - 1]][N - 1];
DET = DET * A[Row[N - 1]][N - 1];

/* Start the forward substitution */
for (K = 0; K < N; K++)
/* Start the forward substitution */
for (K = 0; K < N; K++)
Y[K] = B[K];
Y[0] = B[Row[0]];
for (K = 1; K < N; K++) {
for (K = 1; K < N; K++) {
SUM = 0;
for (C = 0; C <= K - 1; C++)
SUM += A[Row[K]][C] * Y[C];
Y[K] = B[Row[K]] - SUM;
}
}

/* Start the back substitution */
X[N - 1] = Y[N - 1] / A[Row[N - 1]][N - 1];
for (K = N - 2; K >= 0; K--) {
/* Start the back substitution */
X[N - 1] = Y[N - 1] / A[Row[N - 1]][N - 1];
for (K = N - 2; K >= 0; K--) {
SUM = 0;
for (C = K + 1; C < N; C++) {
SUM += A[Row[K]][C] * X[C];
}
X[K] = (Y[K] - SUM) / A[Row[K]][K];
} /* End of back substitution */
} /* End of back substitution */

/* Output */
for (K = 0; K < N; K++)
/* Output */
for (K = 0; K < N; K++)
Cf[K] = X[K];
}

+ 88
- 88
satcal.h View File

@@ -18,107 +18,107 @@
*/

const struct {
float d[4][3];
struct {
float vc, A, B;
} rad[3];
struct {
float Ns;
float b[3];
} cor[3];
float d[4][3];
struct {
float vc, A, B;
} rad[3];
struct {
float Ns;
float b[3];
} cor[3];
/* Calibration corficcients taken from the NOAA KLM satellite user guide
* https://web.archive.org/web/20141220021557/https://www.ncdc.noaa.gov/oa/pod-guide/ncdc/docs/klm/tables.htm
*/
} satcal[] = {
{ // NOAA 15
{ // PRT coefficient d0, d1, d2
{276.60157, 0.051045, 1.36328E-06},
{276.62531, 0.050909, 1.47266E-06},
{276.67413, 0.050907, 1.47656E-06},
{276.59258, 0.050966, 1.47656E-06}
},
{ // Channel radiance coefficient vc, A, B
{925.4075, 0.337810, 0.998719}, // Channel 4
{839.8979, 0.304558, 0.999024}, // Channel 5
{2695.9743, 1.621256, 0.998015} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-4.50, {4.76, -0.0932, 0.0004524}}, // Channel 4
{-3.61, {3.83, -0.0659, 0.0002811}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
{ // PRT coefficient d0, d1, d2
{276.60157, 0.051045, 1.36328E-06},
{276.62531, 0.050909, 1.47266E-06},
{276.67413, 0.050907, 1.47656E-06},
{276.59258, 0.050966, 1.47656E-06}
},
{ // Channel radiance coefficient vc, A, B
{925.4075, 0.337810, 0.998719}, // Channel 4
{839.8979, 0.304558, 0.999024}, // Channel 5
{2695.9743, 1.621256, 0.998015} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-4.50, {4.76, -0.0932, 0.0004524}}, // Channel 4
{-3.61, {3.83, -0.0659, 0.0002811}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
},
{ // NOAA 16
{ // PRT coeff d0, d1, d2
{276.355, 5.562E-02, -1.590E-05},
{276.142, 5.605E-02, -1.707E-05},
{275.996, 5.486E-02, -1.223E-05},
{276.132, 5.494E-02, -1.344E-05}
},
{ // Channel radiance coefficient vc, A, B
{917.2289, 0.332380, 0.998522}, // Channel 4
{838.1255, 0.674623, 0.998363}, // Channel 5
{2700.1148, 1.592459, 0.998147} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-2.467, {2.96, -0.05411, 0.00024532}}, // Channel 4
{-2.009, {2.25, -0.03665, 0.00014854}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
{ // PRT coeff d0, d1, d2
{276.355, 5.562E-02, -1.590E-05},
{276.142, 5.605E-02, -1.707E-05},
{275.996, 5.486E-02, -1.223E-05},
{276.132, 5.494E-02, -1.344E-05}
},
{ // Channel radiance coefficient vc, A, B
{917.2289, 0.332380, 0.998522}, // Channel 4
{838.1255, 0.674623, 0.998363}, // Channel 5
{2700.1148, 1.592459, 0.998147} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-2.467, {2.96, -0.05411, 0.00024532}}, // Channel 4
{-2.009, {2.25, -0.03665, 0.00014854}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
},
{ // NOAA 17
{ // PRT coefficient d0, d1, d2
{276.628, 0.05098, 1.371e-06},
{276.538, 0.05098, 1.371e-06},
{276.761, 0.05097, 1.369e-06},
{276.660, 0.05100, 1.348e-06}
},
{ // Channel radiance coefficient vc, A, B
{926.2947, 0.271683, 0.998794}, // Channel 4
{839.8246, 0.309180, 0.999012}, // Channel 5
{2669.3554, 1.702380, 0.997378} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-8.55, {8.22, -0.15795, 0.00075579}}, // Channel 4
{-3.97, {4.31, -0.07318, 0.00030976}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
{ // PRT coefficient d0, d1, d2
{276.628, 0.05098, 1.371e-06},
{276.538, 0.05098, 1.371e-06},
{276.761, 0.05097, 1.369e-06},
{276.660, 0.05100, 1.348e-06}
},
{ // Channel radiance coefficient vc, A, B
{926.2947, 0.271683, 0.998794}, // Channel 4
{839.8246, 0.309180, 0.999012}, // Channel 5
{2669.3554, 1.702380, 0.997378} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-8.55, {8.22, -0.15795, 0.00075579}}, // Channel 4
{-3.97, {4.31, -0.07318, 0.00030976}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
},
{ // NOAA 18
{ // PRT coefficient d0, d1, d2
{276.601, 0.05090, 1.657e-06},
{276.683, 0.05101, 1.482e-06},
{276.565, 0.05117, 1.313e-06},
{276.615, 0.05103, 1.484e-06}
},
{ // Channel radiance coefficient vc, A, B
{928.1460, 0.436645, 0.998607}, // Channel 4
{833.2532, 0.253179, 0.999057}, // Channel 5
{2659.7952, 1.698704, 0.996960} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-5.53, {5.82, -0.11069, 0.00052337}}, // Channel 4
{-2.22, {2.67, -0.04360, 0.00017715}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
{ // PRT coefficient d0, d1, d2
{276.601, 0.05090, 1.657e-06},
{276.683, 0.05101, 1.482e-06},
{276.565, 0.05117, 1.313e-06},
{276.615, 0.05103, 1.484e-06}
},
{ // Channel radiance coefficient vc, A, B
{928.1460, 0.436645, 0.998607}, // Channel 4
{833.2532, 0.253179, 0.999057}, // Channel 5
{2659.7952, 1.698704, 0.996960} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-5.53, {5.82, -0.11069, 0.00052337}}, // Channel 4
{-2.22, {2.67, -0.04360, 0.00017715}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
},
{ // NOAA 19
{ // PRT coefficient d0, d1, d2
{276.6067, 0.051111, 1.405783E-06},
{276.6119, 0.051090, 1.496037E-06},
{276.6311, 0.051033, 1.496990E-06},
{276.6268, 0.051058, 1.493110E-06}
},
{ // Channel radiance coefficient vc, A, B
{928.9, 0.53959, 0.998534}, // Channel 4
{831.9, 0.36064, 0.998913}, // Channel 5
{2670.0, 1.67396, 0.997364} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-5.49, {5.70 -0.11187, 0.00054668}}, // Channel 4
{-3.39, {3.58 -0.05991, 0.00024985}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
{ // PRT coefficient d0, d1, d2
{276.6067, 0.051111, 1.405783E-06},
{276.6119, 0.051090, 1.496037E-06},
{276.6311, 0.051033, 1.496990E-06},
{276.6268, 0.051058, 1.493110E-06}
},
{ // Channel radiance coefficient vc, A, B
{928.9, 0.53959, 0.998534}, // Channel 4
{831.9, 0.36064, 0.998913}, // Channel 5
{2670.0, 1.67396, 0.997364} // Channel 3B
},
{ // Nonlinear radiance correction Ns, b0, b1, b2
{-5.49, {5.70 -0.11187, 0.00054668}}, // Channel 4
{-3.39, {3.58 -0.05991, 0.00024985}}, // Channel 5
{0.0, {0.0, 0.0, 0.0}} // Channel 3B
}
}};

const float c1 = 1.1910427e-5;

Loading…
Cancel
Save