1 #include <string.h>
2 #include <stdio.h>
3 #include <stdlib.h>
4 #include "fitsio.h"
5
main(int argc,char * argv[])6 int main(int argc, char *argv[])
7 {
8 fitsfile *infptr, *outfptr; /* FITS file pointers defined in fitsio.h */
9 int status = 0, tstatus, ii = 1, iteration = 0, single = 0, hdupos;
10 int hdutype, bitpix, bytepix, naxis = 0, nkeys, datatype = 0, anynul;
11 long naxes[9] = {1, 1, 1, 1, 1, 1, 1, 1, 1};
12 long first, totpix = 0, npix;
13 double *array, bscale = 1.0, bzero = 0.0, nulval = 0.;
14 char card[81];
15
16 if (argc != 3)
17 {
18 printf("\n");
19 printf("Usage: imcopy inputImage outputImage[compress]\n");
20 printf("\n");
21 printf("Copy an input image to an output image, optionally compressing\n");
22 printf("or uncompressing the image in the process. If the [compress]\n");
23 printf("qualifier is appended to the output file name then the input image\n");
24 printf("will be compressed using the tile-compressed format. In this format,\n");
25 printf("the image is divided into rectangular tiles and each tile of pixels\n");
26 printf("is compressed and stored in a variable-length row of a binary table.\n");
27 printf("If the [compress] qualifier is omitted, and the input image is\n");
28 printf("in tile-compressed format, then the output image will be uncompressed.\n");
29 printf("\n");
30 printf("If an extension name or number is appended to the input file name, \n");
31 printf("enclosed in square brackets, then only that single extension will be\n");
32 printf("copied to the output file. Otherwise, every extension in the input file\n");
33 printf("will be processed in turn and copied to the output file.\n");
34 printf("\n");
35 printf("Examples:\n");
36 printf("\n");
37 printf("1) imcopy image.fit 'cimage.fit[compress]'\n");
38 printf("\n");
39 printf(" This compresses the input image using the default parameters, i.e.,\n");
40 printf(" using the Rice compression algorithm and using row by row tiles.\n");
41 printf("\n");
42 printf("2) imcopy cimage.fit image2.fit\n");
43 printf("\n");
44 printf(" This uncompresses the image created in the first example.\n");
45 printf(" image2.fit should be identical to image.fit if the image\n");
46 printf(" has an integer datatype. There will be small differences\n");
47 printf(" in the pixel values if it is a floating point image.\n");
48 printf("\n");
49 printf("3) imcopy image.fit 'cimage.fit[compress GZIP 100,100;q 16]'\n");
50 printf("\n");
51 printf(" This compresses the input image using the following parameters:\n");
52 printf(" GZIP compression algorithm;\n");
53 printf(" 100 X 100 pixel compression tiles;\n");
54 printf(" quantization level = 16 (only used with floating point images)\n");
55 printf("\n");
56 printf("The full syntax of the compression qualifier is:\n");
57 printf(" [compress ALGORITHM TDIM1,TDIM2,...; q QLEVEL s SCALE]\n");
58 printf("where the allowed ALGORITHM values are:\n");
59 printf(" Rice, HCOMPRESS, HSCOMPRESS, GZIP, or PLIO. \n");
60 printf(" (HSCOMPRESS is a variant of HCOMPRESS in which a small\n");
61 printf(" amount of smoothing is applied to the uncompressed image\n");
62 printf(" to help suppress blocky compression artifacts in the image\n");
63 printf(" when using large values for the 'scale' parameter).\n");
64 printf("TDIMn is the size of the compression tile in each dimension,\n");
65 printf("\n");
66 printf("QLEVEL specifies the quantization level when converting a floating\n");
67 printf("point image into integers, prior to compressing the image. The\n");
68 printf("default value = 16, which means the image will be quantized into\n");
69 printf("integer levels that are spaced at intervals of sigma/16., where \n");
70 printf("sigma is the estimated noise level in background areas of the image.\n");
71 printf("If QLEVEL is negative, this means use the absolute value for the\n");
72 printf("quantization spacing (e.g. 'q -0.005' means quantize the floating\n");
73 printf("point image such that the scaled integers represent steps of 0.005\n");
74 printf("in the original image).\n");
75 printf("\n");
76 printf("SCALE is the integer scale factor that only applies to the HCOMPRESS\n");
77 printf("algorithm. The default value SCALE = 0 forces the image to be\n");
78 printf("losslessly compressed; Greater amounts of lossy compression (resulting\n");
79 printf("in smaller compressed files) can be specified with larger SCALE values.\n");
80 printf("\n");
81 printf("\n");
82 printf("Note that it may be necessary to enclose the file names\n");
83 printf("in single quote characters on the Unix command line.\n");
84 return(0);
85 }
86
87 /* Open the input file and create output file */
88 fits_open_file(&infptr, argv[1], READONLY, &status);
89 fits_create_file(&outfptr, argv[2], &status);
90
91 if (status != 0) {
92 fits_report_error(stderr, status);
93 return(status);
94 }
95
96 fits_get_hdu_num(infptr, &hdupos); /* Get the current HDU position */
97
98 /* Copy only a single HDU if a specific extension was given */
99 if (hdupos != 1 || strchr(argv[1], '[')) single = 1;
100
101 for (; !status; hdupos++) /* Main loop through each extension */
102 {
103
104 fits_get_hdu_type(infptr, &hdutype, &status);
105
106 if (hdutype == IMAGE_HDU) {
107
108 /* get image dimensions and total number of pixels in image */
109 for (ii = 0; ii < 9; ii++)
110 naxes[ii] = 1;
111
112 fits_get_img_param(infptr, 9, &bitpix, &naxis, naxes, &status);
113
114 totpix = naxes[0] * naxes[1] * naxes[2] * naxes[3] * naxes[4]
115 * naxes[5] * naxes[6] * naxes[7] * naxes[8];
116 }
117
118 if (hdutype != IMAGE_HDU || naxis == 0 || totpix == 0) {
119
120 /* just copy tables and null images */
121 fits_copy_hdu(infptr, outfptr, 0, &status);
122
123 } else {
124
125 /* Explicitly create new image, to support compression */
126 fits_create_img(outfptr, bitpix, naxis, naxes, &status);
127 if (status) {
128 fits_report_error(stderr, status);
129 return(status);
130 }
131
132 if (fits_is_compressed_image(outfptr, &status)) {
133 /* write default EXTNAME keyword if it doesn't already exist */
134 tstatus = 0;
135 fits_read_card(infptr, "EXTNAME", card, &tstatus);
136 if (tstatus) {
137 strcpy(card, "EXTNAME = 'COMPRESSED_IMAGE' / name of this binary table extension");
138 fits_write_record(outfptr, card, &status);
139 }
140 }
141
142 /* copy all the user keywords (not the structural keywords) */
143 fits_get_hdrspace(infptr, &nkeys, NULL, &status);
144
145 for (ii = 1; ii <= nkeys; ii++) {
146 fits_read_record(infptr, ii, card, &status);
147 if (fits_get_keyclass(card) > TYP_CMPRS_KEY)
148 fits_write_record(outfptr, card, &status);
149 }
150
151 /* delete default EXTNAME keyword if it exists */
152 /*
153 if (!fits_is_compressed_image(outfptr, &status)) {
154 tstatus = 0;
155 fits_read_key(outfptr, TSTRING, "EXTNAME", card, NULL, &tstatus);
156 if (!tstatus) {
157 if (strcmp(card, "COMPRESSED_IMAGE") == 0)
158 fits_delete_key(outfptr, "EXTNAME", &status);
159 }
160 }
161 */
162
163 switch(bitpix) {
164 case BYTE_IMG:
165 datatype = TBYTE;
166 break;
167 case SHORT_IMG:
168 datatype = TSHORT;
169 break;
170 case LONG_IMG:
171 datatype = TINT;
172 break;
173 case FLOAT_IMG:
174 datatype = TFLOAT;
175 break;
176 case DOUBLE_IMG:
177 datatype = TDOUBLE;
178 break;
179 }
180
181 bytepix = abs(bitpix) / 8;
182
183 npix = totpix;
184 iteration = 0;
185
186 /* try to allocate memory for the entire image */
187 /* use double type to force memory alignment */
188 array = (double *) calloc(npix, bytepix);
189
190 /* if allocation failed, divide size by 2 and try again */
191 while (!array && iteration < 10) {
192 iteration++;
193 npix = npix / 2;
194 array = (double *) calloc(npix, bytepix);
195 }
196
197 if (!array) {
198 printf("Memory allocation error\n");
199 return(0);
200 }
201
202 /* turn off any scaling so that we copy the raw pixel values */
203 fits_set_bscale(infptr, bscale, bzero, &status);
204 fits_set_bscale(outfptr, bscale, bzero, &status);
205
206 first = 1;
207 while (totpix > 0 && !status)
208 {
209 /* read all or part of image then write it back to the output file */
210 fits_read_img(infptr, datatype, first, npix,
211 &nulval, array, &anynul, &status);
212
213 fits_write_img(outfptr, datatype, first, npix, array, &status);
214 totpix = totpix - npix;
215 first = first + npix;
216 }
217 free(array);
218 }
219
220 if (single) break; /* quit if only copying a single HDU */
221 fits_movrel_hdu(infptr, 1, NULL, &status); /* try to move to next HDU */
222 }
223
224 if (status == END_OF_FILE) status = 0; /* Reset after normal error */
225
226 fits_close_file(outfptr, &status);
227 fits_close_file(infptr, &status);
228
229 /* if error occurred, print out error message */
230 if (status)
231 fits_report_error(stderr, status);
232 return(status);
233 }
234