You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

main.c 12 KiB

21 years ago
21 years ago
21 years ago
21 years ago
21 years ago
19 years ago
21 years ago
21 years ago
19 years ago
21 years ago
21 years ago
21 years ago
21 years ago
21 years ago
21 years ago
19 years ago
21 years ago
21 years ago
21 years ago
21 years ago
21 years ago
21 years ago
21 years ago
123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446
  1. /*
  2. * This file is part of Aptdec.
  3. * Copyright (c) 2004-2009 Thierry Leconte (F4DWV), Xerbo (xerbo@protonmail.com) 2019
  4. *
  5. * Aptdec is free software: you can redistribute it and/or modify
  6. * it under the terms of the GNU General Public License as published by
  7. * the Free Software Foundation, either version 2 of the License, or
  8. * (at your option) any later version.
  9. *
  10. * This program is distributed in the hope that it will be useful,
  11. * but WITHOUT ANY WARRANTY; without even the implied warranty of
  12. * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  13. * GNU General Public License for more details.
  14. *
  15. * You should have received a copy of the GNU General Public License
  16. * along with this program. If not, see <https://www.gnu.org/licenses/>.
  17. *
  18. */
  19. #include <stdlib.h>
  20. #include <stdio.h>
  21. #include <string.h>
  22. #include <getopt.h>
  23. #include <libgen.h>
  24. #include <math.h>
  25. #include <sndfile.h>
  26. #include <png.h>
  27. #include "messages.h"
  28. #include "offsets.h"
  29. #include "palette.h"
  30. extern int getpixelrow(float *pixelv);
  31. extern int init_dsp(double F);
  32. static SNDFILE *inwav;
  33. static int initsnd(char *filename) {
  34. SF_INFO infwav;
  35. int res;
  36. // Open audio file
  37. infwav.format = 0;
  38. inwav = sf_open(filename, SFM_READ, &infwav);
  39. if (inwav == NULL) {
  40. fprintf(stderr, ERR_FILE_READ, filename);
  41. return(0);
  42. }
  43. res = init_dsp(infwav.samplerate);
  44. printf("Input file: %s\n", filename);
  45. if(res < 0) {
  46. fprintf(stderr, "Input sample rate too low: %d\n", infwav.samplerate);
  47. return(0);
  48. }else if(res > 0) {
  49. fprintf(stderr, "Input sample rate too high: %d\n", infwav.samplerate);
  50. return(0);
  51. }
  52. printf("Input sample rate: %d\n", infwav.samplerate);
  53. if (infwav.channels != 1) {
  54. fprintf(stderr, "Too many channels in input file: %d\n", infwav.channels);
  55. return(0);
  56. }
  57. return(1);
  58. }
  59. // Get samples from the wave file
  60. int getsample(float *sample, int nb) {
  61. return sf_read_float(inwav, sample, nb);
  62. }
  63. static png_text text_ptr[] = {
  64. {PNG_TEXT_COMPRESSION_NONE, "Software", VERSION},
  65. {PNG_TEXT_COMPRESSION_NONE, "Channel", NULL, 0},
  66. {PNG_TEXT_COMPRESSION_NONE, "Description", "NOAA satellite image", 20}
  67. };
  68. // TODO: this function needs to be tidied up
  69. /* Effects
  70. * 0 - Nothing
  71. * 1 - Crop telemetry
  72. * 2 - False color
  73. * 3 - Layered
  74. */
  75. static int ImageOut(char *filename, char *chid, float **prow, int nrow, int width, int offset, png_color *palette, int effects) {
  76. FILE *pngfile;
  77. png_infop info_ptr;
  78. png_structp png_ptr;
  79. // Reduce the width of the image to componsate for the missing telemetry
  80. if(effects == 1) width -= TOTAL_TELE;
  81. // Initalise the PNG writer
  82. png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);
  83. if (!png_ptr) {
  84. fprintf(stderr, ERR_PNG_WRITE);
  85. return(0);
  86. }
  87. // Metadata
  88. info_ptr = png_create_info_struct(png_ptr);
  89. if (!info_ptr) {
  90. png_destroy_write_struct(&png_ptr, (png_infopp) NULL);
  91. fprintf(stderr, ERR_PNG_INFO);
  92. return(0);
  93. }
  94. extern void falsecolor(float vis, float temp, float *r, float *g, float *b);
  95. if(effects == 2){
  96. // 8 bit RGB image
  97. png_set_IHDR(png_ptr, info_ptr, width, nrow,
  98. 8, PNG_COLOR_TYPE_RGB, PNG_INTERLACE_NONE,
  99. PNG_COMPRESSION_TYPE_DEFAULT, PNG_FILTER_TYPE_DEFAULT);
  100. }else if(palette == NULL) {
  101. // Greyscale image
  102. png_set_IHDR(png_ptr, info_ptr, width, nrow,
  103. 8, PNG_COLOR_TYPE_GRAY, PNG_INTERLACE_NONE,
  104. PNG_COMPRESSION_TYPE_DEFAULT, PNG_FILTER_TYPE_DEFAULT);
  105. } else {
  106. // Palleted color image
  107. png_set_IHDR(png_ptr, info_ptr, width, nrow,
  108. 8, PNG_COLOR_TYPE_PALETTE, PNG_INTERLACE_NONE,
  109. PNG_COMPRESSION_TYPE_DEFAULT, PNG_FILTER_TYPE_DEFAULT);
  110. png_set_PLTE(png_ptr, info_ptr, palette, 256);
  111. }
  112. text_ptr[1].text = chid;
  113. text_ptr[1].text_length = strlen(chid);
  114. png_set_text(png_ptr, info_ptr, text_ptr, 3);
  115. png_set_pHYs(png_ptr, info_ptr, 4000, 4000, PNG_RESOLUTION_METER);
  116. if(effects == 2){
  117. printf("Computing false color & writing: %s", filename);
  118. }else{
  119. printf("Writing %s", filename);
  120. }
  121. fflush(stdout);
  122. pngfile = fopen(filename, "wb");
  123. if (pngfile == NULL) {
  124. fprintf(stderr, ERR_FILE_WRITE, filename);
  125. return(1);
  126. }
  127. png_init_io(png_ptr, pngfile);
  128. png_write_info(png_ptr, info_ptr);
  129. for (int n = 0; n < nrow; n++) {
  130. png_color pix[width];
  131. png_byte pixel[width];
  132. float *pixelv = prow[n];
  133. int f = 0;
  134. for (int i = 0; i < width; i++) {
  135. // Skip parts of the image that are telemetry
  136. if(effects == 1){
  137. switch (i) {
  138. case 0:
  139. f += SYNC_WIDTH + SPC_WIDTH;
  140. break;
  141. case CH_WIDTH:
  142. f += TELE_WIDTH + SYNC_WIDTH + SPC_WIDTH;
  143. break;
  144. case CH_WIDTH*2:
  145. f += TELE_WIDTH;
  146. }
  147. }
  148. if(effects == 2){
  149. float r = 0, g = 0, b = 0;
  150. falsecolor(pixelv[i+CHA_OFFSET], pixelv[i+CHB_OFFSET], &r, &g, &b);
  151. pix[i].red = r;
  152. pix[i].green = g;
  153. pix[i].blue = b;
  154. }else if(effects == 3){
  155. // Layered image, overlay clouds in channel B over channel A
  156. float cloud = CLIP(pixelv[i+CHB_OFFSET]-141, 0, 255)/114;
  157. pixel[i] = MCOMPOSITE(240, cloud, pixelv[i+CHA_OFFSET], 1);
  158. }else{
  159. pixel[i] = CLIP(pixelv[i + offset + f], 0, 255);
  160. }
  161. }
  162. if(effects == 2){
  163. png_write_row(png_ptr, (png_bytep) pix);
  164. }else{
  165. png_write_row(png_ptr, pixel);
  166. }
  167. }
  168. png_write_end(png_ptr, info_ptr);
  169. fclose(pngfile);
  170. printf("\nDone\n");
  171. png_destroy_write_struct(&png_ptr, &info_ptr);
  172. return(1);
  173. }
  174. // Outputs a image with the value distribution between channel A and B
  175. static void distrib(char *filename, float **prow, int nrow) {
  176. float *distrib[256];
  177. int max = 0;
  178. // Assign memory
  179. for(int i = 0; i < 256; i++)
  180. distrib[i] = (float *) malloc(sizeof(float) * 256);
  181. for(int n = 0; n < nrow; n++) {
  182. float *pixelv = prow[n];
  183. for(int i = 0; i < CH_WIDTH; i++) {
  184. int y = (int)(pixelv[i + CHA_OFFSET]);
  185. int x = (int)(pixelv[i + CHB_OFFSET]);
  186. distrib[y][x] += 1;
  187. if(distrib[y][x] > max) max = distrib[y][x];
  188. }
  189. }
  190. // Scale to 0-255
  191. for(int x = 0; x < 256; x++)
  192. for(int y = 0; y < 256; y++)
  193. distrib[y][x] = distrib[y][x] / max * 255;
  194. ImageOut(filename, "Brightness distribution", distrib, 256, 256, 0, NULL, 0);
  195. }
  196. extern int calibrate(float **prow, int nrow, int offset, int width, int calibrate);
  197. extern void histogramEqualise(float **prow, int nrow, int offset, int width);
  198. extern void temperature(float **prow, int nrow, int ch, int offset);
  199. extern int Ngvi(float **prow, int nrow);
  200. extern void readfcconf(char *file);
  201. extern int optind;
  202. extern char *optarg;
  203. // Default to NOAA 19
  204. int satnum = 4;
  205. static void usage(void) {
  206. printf("Aptdec [options] audio files ...\n"
  207. "Options:\n"
  208. " -e [c|t] Enhancements\n"
  209. " c: Contrast calibration\n"
  210. " t: Crop telemetry\n"
  211. " h: Histogram equalise\n"
  212. " -i [r|a|b|c|t] Output image type\n"
  213. " r: Raw\n"
  214. " a: Channel A\n"
  215. " b: Channel B\n"
  216. " c: False color\n"
  217. " t: Temperature\n"
  218. " l: Layered\n"
  219. " -d <dir> Image destination directory.\n"
  220. " -s [15-19] Satellite number\n"
  221. " -c <file> False color config file\n");
  222. exit(1);
  223. }
  224. int main(int argc, char **argv) {
  225. char pngfilename[1024];
  226. char name[500];
  227. char pngdirname[500] = "";
  228. // Default to a raw image, with equalization and cropped telemetry
  229. char imgopt[20] = "r";
  230. char enchancements[20] = "ct";
  231. // Image buffer
  232. float *prow[3000];
  233. int nrow;
  234. // Mapping between telemetry wedge value and channel
  235. static struct {
  236. char *id[7];
  237. char *name[7];
  238. } ch = {
  239. { "?", "1", "2", "3A", "4", "5", "3B" },
  240. { "unknown", "visble", "near-infrared", "mid-infrared", "thermal-infrared", "thermal-infrared", "mid-infrared" }
  241. };
  242. // The active sensor in each channel
  243. int chA, chB;
  244. // Print version
  245. printf(VERSION"\n");
  246. // Print usage if there are no arguments
  247. if(argc == 1)
  248. usage();
  249. int c;
  250. while ((c = getopt(argc, argv, "c:d:i:s:e:")) != EOF) {
  251. switch (c) {
  252. // Output directory name
  253. case 'd':
  254. strcpy(pngdirname, optarg);
  255. break;
  256. // False color config file
  257. case 'c':
  258. readfcconf(optarg);
  259. break;
  260. // Output image type
  261. case 'i':
  262. strncpy(imgopt, optarg, 20);
  263. break;
  264. // Satellite number (for calibration)
  265. case 's':
  266. satnum = atoi(optarg)-15;
  267. // Check if it's within the valid range
  268. if (satnum < 0 || satnum > 4) {
  269. fprintf(stderr, "Invalid satellite number, it must be the range [15-19]\n");
  270. exit(1);
  271. }
  272. break;
  273. // Enchancements
  274. case 'e':
  275. strncpy(enchancements, optarg, 20);
  276. break;
  277. default:
  278. usage();
  279. }
  280. }
  281. if(optind == argc){
  282. printf("No audio files provided.\n");
  283. usage();
  284. }
  285. // Process the provided files
  286. for (; optind < argc; optind++) {
  287. chA = chB = 0;
  288. // Generate output name
  289. strcpy(pngfilename, argv[optind]);
  290. strcpy(name, basename(pngfilename));
  291. strtok(name, ".");
  292. if (pngdirname[0] == '\0')
  293. strcpy(pngdirname, dirname(pngfilename));
  294. // Open sound file, exit if that fails
  295. if (initsnd(argv[optind]) == 0)
  296. exit(1);
  297. // Main image building loop
  298. for (nrow = 0; nrow < 3000; nrow++) {
  299. // Allocate 2150 floats worth of memory for every line of the image
  300. prow[nrow] = (float *) malloc(sizeof(float) * 2150);
  301. // Read into prow and break the loop once we reach the end of the image
  302. if (getpixelrow(prow[nrow]) == 0)
  303. break;
  304. printf("Row: %d\r", nrow);
  305. fflush(stdout);
  306. }
  307. printf("\nTotal rows: %d\n", nrow);
  308. sf_close(inwav);
  309. chA = calibrate(prow, nrow, CHA_OFFSET, CH_WIDTH, 0);
  310. chB = calibrate(prow, nrow, CHB_OFFSET, CH_WIDTH, 0);
  311. printf("Channel A: %s (%s)\n", ch.id[chA], ch.name[chA]);
  312. printf("Channel B: %s (%s)\n", ch.id[chB], ch.name[chB]);
  313. // Temperature
  314. if (CONTAINS(imgopt, 't') && chB >= 4) {
  315. temperature(prow, nrow, chB, CHB_OFFSET);
  316. sprintf(pngfilename, "%s/%s-t.png", pngdirname, name);
  317. ImageOut(pngfilename, "Temperature", prow, nrow, CH_WIDTH, CHB_OFFSET, (png_color*)TempPalette, 0);
  318. }
  319. // Run the brightness calibration here because the temperature calibration requires raw data
  320. // Layered & false color images both also need brightness calibration
  321. if(CONTAINS(enchancements, 'c') || CONTAINS(enchancements, 'h') || CONTAINS(imgopt, 'l') || CONTAINS(imgopt, 'c'))
  322. calibrate(prow, nrow, CHA_OFFSET, CH_WIDTH+TELE_WIDTH+SYNC_WIDTH+SPC_WIDTH+CH_WIDTH, 1);
  323. // Histogram equalise
  324. if(CONTAINS(enchancements, 'h')){
  325. histogramEqualise(prow, nrow, CHA_OFFSET, CH_WIDTH);
  326. histogramEqualise(prow, nrow, CHB_OFFSET, CH_WIDTH);
  327. }
  328. // Layered
  329. if (CONTAINS(imgopt, 'l')){
  330. if(chA == 1){
  331. sprintf(pngfilename, "%s/%s-l.png", pngdirname, name);
  332. ImageOut(pngfilename, "Layered", prow, nrow, CH_WIDTH, 0, NULL, 3);
  333. }else{
  334. fprintf(stderr, "Lacking channels required for generting a layered image.\n");
  335. }
  336. }
  337. // Raw image
  338. if (CONTAINS(imgopt, 'r')) {
  339. int croptele = CONTAINS(enchancements, 't');
  340. char channelstr[45];
  341. sprintf(channelstr, "%s (%s) & %s (%s)", ch.id[chA], ch.name[chA], ch.id[chB], ch.name[chB]);
  342. sprintf(pngfilename, "%s/%s-r.png", pngdirname, name);
  343. ImageOut(pngfilename, channelstr, prow, nrow, IMG_WIDTH, 0, NULL, croptele ? 1 : 0);
  344. }
  345. // Channel A
  346. if (CONTAINS(imgopt, 'a')) {
  347. char channelstr[21];
  348. sprintf(channelstr, "%s (%s)", ch.id[chA], ch.name[chA]);
  349. sprintf(pngfilename, "%s/%s-%s.png", pngdirname, name, ch.id[chA]);
  350. ImageOut(pngfilename, channelstr, prow, nrow, CH_WIDTH, CHA_OFFSET, NULL, 0);
  351. }
  352. // Channel B
  353. if (CONTAINS(imgopt, 'b')) {
  354. char channelstr[21];
  355. sprintf(channelstr, "%s (%s)", ch.id[chB], ch.name[chB]);
  356. sprintf(pngfilename, "%s/%s-%s.png", pngdirname, name, ch.id[chB]);
  357. ImageOut(pngfilename, channelstr, prow, nrow, CH_WIDTH , CHB_OFFSET, NULL, 0);
  358. }
  359. // Distribution image
  360. if (CONTAINS(imgopt, 'd')) {
  361. sprintf(pngfilename, "%s/%s-d.png", pngdirname, name);
  362. distrib(pngfilename, prow, nrow);
  363. }
  364. // False color image
  365. if(CONTAINS(imgopt, 'c')){
  366. if (chA == 2 && chB >= 4) { // Normal false color
  367. sprintf(pngfilename, "%s/%s-c.png", pngdirname, name);
  368. //ImageRGBOut(pngfilename, prow, nrow);
  369. ImageOut(pngfilename, "False Color", prow, nrow, CH_WIDTH, 0, NULL, 2);
  370. } else if (chB == 2) { // GVI (global vegetation index) false color
  371. Ngvi(prow, nrow);
  372. sprintf(pngfilename, "%s/%s-c.png", pngdirname, name);
  373. ImageOut(pngfilename, "GVI False Color", prow, nrow, CH_WIDTH, CHB_OFFSET, (png_color*)GviPalette, 0);
  374. } else {
  375. fprintf(stderr, "Lacking channels required for generating a false color image.\n");
  376. }
  377. }
  378. }
  379. exit(0);
  380. }