add save to FITS

Makefile

@@ -1,6 +1,7 @@
 # run `make DEF=...` to add extra defines
 PROGRAM := sbig340
-LDFLAGS := -ltiff -lm -fdata-sections -ffunction-sections -Wl,--gc-sections -Wl,--discard-all
+LDFLAGS := -fdata-sections -ffunction-sections -Wl,--gc-sections -Wl,--discard-all
+LDFLAGS += -ltiff -lm $(shell pkg-config --libs cfitsio)
 SRCS := $(wildcard *.c)
 DEFINES := $(DEF) -D_XOPEN_SOURCE=1111 -DEBUG
 OBJDIR := mk

imfunctions.c

@@ -26,7 +26,10 @@
 #include "term.h"
 #include <strings.h> // strncasecmp
 #include <tiffio.h>  // save tiff
-#include <math.h>
+#include <math.h>    // sqrt
+#include <time.h>    // time, gmtime etc
+#include <fitsio.h>  // save fits
+#include <libgen.h>  // basename
 
 // image type suffixes
 #define SUFFIX_FITS         "fits"
@@ -50,7 +53,11 @@ static char *make_filename(const char *outfile, const char *suff, store_type st)
             }
             else return buff;
         }else{ // file exists
-            if(st == STORE_REWRITE) return buff;
+            if(st == STORE_REWRITE){
+                if(0 == strcmp(suff, SUFFIX_FITS)) // add '!' before image name
+                snprintf(buff, FILENAME_MAX, "!%s.%s", outfile, suff);
+                return buff;
+            }
             else return NULL; // file exists on option STORE_NORMAL
         }
     }
@@ -87,33 +94,34 @@ imstorage *chk_storeimg(char *filename, char* store, char *format){
     char *nm = strdup(filename);
     if(!nm) ERRX("strdup");
     char *pt = strrchr(nm, '.');
+    image_format fbysuff = FORMAT_NONE;
+    // check if name's suffix is filetype
     if(pt){
         char *suff = pt + 1;
         DBG("got suffix: %s", suff);
-        // check if name's suffix is filetype
-        image_format fbysuff = FORMAT_NONE;
         if(0 == strcasecmp(suff, "tiff") || 0 == strcasecmp(suff, "tif")) fbysuff = FORMAT_TIFF;
         else if(0 == strcasecmp(suff, "fits") || 0 == strcasecmp(suff, "fit")) fbysuff = FORMAT_FITS;
         else if(0 == strcasecmp(suff, "raw") || 0 == strcasecmp(suff, "bin") || 0 == strcasecmp(suff, "dump"))
             fbysuff = FORMAT_RAW;
         if(fbysuff != FORMAT_NONE) *pt = 0; // truncate filename if suffix found
-        if(format){ // check if user gave format
-            fbysuff = FORMAT_NONE;
-            if(strchr(format, 'f') || strchr(format, 'F')) fbysuff |= FORMAT_FITS;
-            if(strchr(format, 't') || strchr(format, 'T')) fbysuff |= FORMAT_TIFF;
-            if(strchr(format, 'r') || strchr(format, 'R') ||
-               strchr(format, 'd') || strchr(format, 'D')) fbysuff |= FORMAT_RAW;
-            if(fbysuff == FORMAT_NONE){
-                WARNX(_("Wrong format string: %s"), format);
-                free(nm);
-                return NULL;
-            }
-            fmt = fbysuff;
-        }else{// if user gave image with suffix, change format; else leave default: FITS
-            if(fbysuff != FORMAT_NONE) fmt = fbysuff;
-            DBG("fmt: %d", fmt);
-        }
     }
+    if(format){ // check if user gave format
+        fbysuff = FORMAT_NONE;
+        if(strchr(format, 'f') || strchr(format, 'F')) fbysuff |= FORMAT_FITS;
+        if(strchr(format, 't') || strchr(format, 'T')) fbysuff |= FORMAT_TIFF;
+        if(strchr(format, 'r') || strchr(format, 'R') ||
+           strchr(format, 'd') || strchr(format, 'D')) fbysuff |= FORMAT_RAW;
+        if(fbysuff == FORMAT_NONE){
+            WARNX(_("Wrong format string: %s"), format);
+            free(nm);
+            return NULL;
+        }
+        fmt = fbysuff;
+    }else{// if user gave image with suffix, change format; else leave default: FITS
+        if(fbysuff != FORMAT_NONE) fmt = fbysuff;
+        DBG("fmt: %d", fmt);
+    }
+
     // now check all names
     #define FMTSZ (3)
     image_format formats[FMTSZ] = {FORMAT_FITS, FORMAT_TIFF, FORMAT_RAW};
@@ -135,7 +143,7 @@ imstorage *chk_storeimg(char *filename, char* store, char *format){
 
 /**
  * Try to write tiff file
- * @return 1 if all OK
+ * @return 0 if all OK
  */
 int writetiff(imstorage *img){
     int H = img->H, W = img->W, y;
@@ -144,7 +152,7 @@ int writetiff(imstorage *img){
     TIFF *image = NULL;
     if(!name || !(image = TIFFOpen(name, "w"))){
         WARN("Can't save tiff file");
-        return 0;
+        return 1;
     }
     TIFFSetField(image, TIFFTAG_IMAGEWIDTH, W);
     TIFFSetField(image, TIFFTAG_IMAGELENGTH, H);
@@ -161,21 +169,25 @@ int writetiff(imstorage *img){
         if(S != (size_t)TIFFWriteEncodedStrip(image, y, data, S)){
             WARNX(_("Error writing %s"), name);
             TIFFClose(image);
-            return 0;
+            return 2;
         }
         //TIFFWriteScanline(image, data, y, 0);
     }
     TIFFClose(image);
     green(_("Image %s saved\n"), name);
-    return 1;
+    return 0;
 }
 
+/**
+ * Calculate image statistics: print it on screen and save for `writefits`
+ */
+static uint16_t glob_min, glob_max, glob_avr, glob_std;
 void print_stat(imstorage *img){
     size_t size = img->W*img->H, i, Noverld = 0L, N = 0L;
     double pv, sum=0., sum2=0., sz = (double)size, tres;
     uint16_t *ptr = img->imdata, val, valoverld;
     uint16_t max = 0, min = 65535;
-    valoverld = min - 5;
+    valoverld = 65530;
     for(i = 0; i < size; i++, ptr++){
         val = *ptr;
         pv = (double) val;
@@ -187,6 +199,7 @@ void print_stat(imstorage *img){
     }
     printf(_("Image stat:\n"));
     double avr = sum/sz, std = sqrt(fabs(sum2/sz - avr*avr));
+    glob_avr = avr, glob_std = std, glob_max = max, glob_min = min;
     printf("avr = %.1f, std = %.1f, Noverload = %ld\n", avr, std, Noverld);
     printf("max = %u, min = %u, W*H = %ld\n", max, min, size);
     Noverld = 0L;
@@ -212,6 +225,98 @@ void print_stat(imstorage *img){
     printf("At 3sigma: Noverload = %ld, avr = %.3f, std = %.3f\n", Noverld, avr, std);
 }
 
+#define TRYFITS(f, ...)                     \
+do{ int status = 0;                         \
+    f(__VA_ARGS__, &status);                \
+    if (status){                            \
+        fits_report_error(stderr, status);  \
+        return -1;}                         \
+}while(0)
+#define WRITEKEY(...)                           \
+do{ int status = 0;                             \
+    fits_write_key(fp, __VA_ARGS__, &status);       \
+    if(status) fits_report_error(stderr, status);\
+}while(0)
+
+int writefits(imstorage *img){
+    long naxes[2] = {img->W, img->H};
+    char buf[80];
+    fitsfile *fp;
+    char *filename = make_filename(img->imname, SUFFIX_FITS, img->st);
+    TRYFITS(fits_create_file, &fp, filename);
+    TRYFITS(fits_create_img, fp, USHORT_IMG, 2, naxes);
+    // FILE / Input file original name
+    char *fn = basename(filename);
+    if(*fn == '!') ++fn; // remove '!' from filename
+    WRITEKEY(TSTRING, "FILE", fn, "Input file original name");
+    // ORIGIN / organization responsible for the data
+    WRITEKEY(TSTRING, "ORIGIN", "SAO RAS", "organization responsible for the data");
+    // OBSERVAT / Observatory name
+    WRITEKEY(TSTRING, "OBSERVAT", "Special Astrophysical Observatory, Russia", "Observatory name");
+    // DETECTOR / detector
+    WRITEKEY(TSTRING, "DETECTOR", "Kodak KAI-340", "Detector model");
+    // INSTRUME / Instrument
+    WRITEKEY(TSTRING, "INSTRUME", "SBIG All-sky 340C", "Instrument");
+    // PXSIZE / pixel size
+    WRITEKEY(TSTRING, "PXSIZE", "7.4 x 7.4", "Pixel size in um");
+    WRITEKEY(TSTRING, "FIELD", "180 degrees", "Camera field of view");
+    switch(img->imtype){
+        case IMTYPE_AUTODARK:
+            sprintf(buf, "object without dark");
+        break;
+        case IMTYPE_DARK:
+            sprintf(buf, "dark");
+        break;
+        case IMTYPE_LIGHT:
+        default:
+            sprintf(buf, "object");
+    }
+    // IMAGETYP / object, flat, dark, bias, scan, eta, neon, push
+    WRITEKEY(TSTRING, "IMAGETYP", buf, "Image type");
+    // DATAMAX, DATAMIN / Max,min pixel value
+    uint16_t val = UINT16_MAX;
+    WRITEKEY(TUSHORT, "DATAMAX", &val, "Max pixel value");
+    val = 0;
+    WRITEKEY(TUSHORT, "DATAMIN", &val, "Min pixel value");
+    // statistical values
+    WRITEKEY(TUSHORT, "STATMAX", &glob_max, "Max pixel value");
+    WRITEKEY(TUSHORT, "STATMIN", &glob_min, "Min pixel value");
+    WRITEKEY(TUSHORT, "STATAVR", &glob_avr, "Average pixel value");
+    WRITEKEY(TUSHORT, "STATSTD", &glob_std, "Standart deviation of pixel value");
+    // EXPTIME / actual exposition time (sec)
+    WRITEKEY(TDOUBLE, "EXPTIME", &img->exptime, "actual exposition time (sec)");
+    // DATE / Creation date (YYYY-MM-DDThh:mm:ss, UTC)
+    time_t savetime = time(NULL);
+    strftime(buf, 79, "%Y-%m-%dT%H:%M:%S", gmtime(&savetime));
+    WRITEKEY(TSTRING, "DATE", buf, "Creation date (YYYY-MM-DDThh:mm:ss, UTC)");
+    struct tm *tm_starttime = gmtime(&img->exposetime);
+    strftime(buf, 79, "exposition starts at %d/%m/%Y, %H:%M:%S (UTC)", tm_starttime);
+    //long tstart = (long)img->exposetime;
+    WRITEKEY(TLONG, "UNIXTIME", &img->exposetime, buf);
+    tm_starttime = localtime(&img->exposetime);
+    strftime(buf, 79, "%Y/%m/%d", tm_starttime);
+    // DATE-OBS / DATE (YYYY/MM/DD) OF OBS.
+    WRITEKEY(TSTRING, "DATE-OBS", buf, "DATE OF OBS. (YYYY/MM/DD, local)");
+    strftime(buf, 79, "%H:%M:%S", tm_starttime);
+    // START / Measurement start time (local) (hh:mm:ss)
+    WRITEKEY(TSTRING, "START", buf, "Measurement start time (hh:mm:ss, local)");
+    // OBJECT  / Object name
+    WRITEKEY(TSTRING, "OBJECT", "sky", "Object name");
+    // BINNING / Binning
+    if(img->binning == 2){
+        WRITEKEY(TSTRING, "BINNING", "2 x 2", "Binning (hbin x vbin)");
+    }
+    if(img->subframe){ // subframe
+        snprintf(buf, 80, "(%d, %d)", img->subframe->Xstart, img->subframe->Ystart);
+        WRITEKEY(TSTRING, "SUBFRAME", buf, "Subframe start coordinates (Xstart, Ystart)");
+    }
+    TRYFITS(fits_write_img, fp, TUSHORT, 1, img->W * img->H, img->imdata);
+    TRYFITS(fits_close_file, fp);
+    if(*filename == '!') ++filename; // remove '!' from filename
+    green(_("Image %s saved\n"), filename);
+    return 0;
+}
+
 
 /**
  * Receive image data & fill img->imdata
@@ -234,17 +339,19 @@ uint16_t *get_imdata(imstorage *img){
  * @return 0 if all OK
  */
 int save_histo(FILE *f, imstorage *img){
-    uint8_t histogram[256];
     if(!img->imdata) return 1000;
+    uint8_t histogram[256];
     size_t l, S = img->W*img->H;
     uint16_t *ptr = img->imdata;
+    memset(histogram, 0, 256);
     for(l = 0; l < S; ++l, ++ptr){
-        ++histogram[(*ptr>>8)&0xff];
+        ++histogram[((*ptr)>>8)&0xff];
     }
     for(l = 0; l < 256; ++l){
         int status = fprintf(f, "%zd\t%u\n", l, histogram[l]);
-        if(status < 0)
+        if(status < 0){
             return status;
+        }
     }
     return 0;
 }
@@ -265,7 +372,7 @@ int writedump(imstorage *img){
         WARN(_("Can't make dump"));
         return 3;
     }
-    name = make_filename(img->imname, "_histogram.txt", img->st);
+    name = make_filename(img->imname, "histogram.txt", img->st);
     if(!name) return 4;
     FILE *h = fopen(name, "w");
     if(!h) return 5;
@@ -290,15 +397,16 @@ int writedump(imstorage *img){
  */
 int store_image(imstorage *img){
     int status = 0;
+    print_stat(img);
     if(!img->imdata && !get_imdata(img)) return 1;
     if(img->imformat & FORMAT_TIFF){ // save tiff file
-        if(!writetiff(img)) status |= 1;
+        if(writetiff(img)) status |= 1;
     }
     if(img->imformat & FORMAT_RAW){
-        if(!writedump(img)) status |= 2;
+        if(writedump(img)) status |= 2;
     }
     if(img->imformat & FORMAT_FITS){ // not supported yet
-        status |= 4;
+        if(writefits(img)) status |= 4;
     }
     return status;
 }

imfunctions.h

@@ -63,8 +63,9 @@ typedef struct{
     double exptime;
     int binning;
     imsubframe *subframe;
-    size_t W, H;      // image size
+    size_t W, H;       // image size
-    uint16_t *imdata; // image data itself
+    uint16_t *imdata;  // image data itself
+    time_t exposetime; // time of exposition start
 } imstorage;
 
 imstorage *chk_storeimg(char *filename, char* store, char *format);

main.c

@@ -75,7 +75,6 @@ int main(int argc, char **argv){
                 if(!get_imdata(img)){
                     WARNX(_("Error image transfer"));
                 }else{
-                    print_stat(img);
                     if(store_image(img))
                         WARNX(_("Error storing image"));
                 }

term.c

@@ -21,6 +21,7 @@
 #include "usefull_macros.h"
 #include "term.h"
 #include <strings.h> // strncasecmp
+#include <time.h>    // time(NULL)
 
 #define BUFLEN 1024
 
@@ -543,6 +544,7 @@ int start_exposition(imstorage *im, char *imtype){
             H = IMHEIGHT;
     }
     im->W = W; im->H = H;
+    im->exposetime = time(NULL);
     return 0;
 }