1 /* Note: this file has been modified from its original form by the
2 Astrometry.net team. For details see http://astrometry.net */
3
4 /* $Id: qfits_image.c,v 1.12 2006/02/23 11:25:25 yjung Exp $
5 *
6 * This file is part of the ESO QFITS Library
7 * Copyright (C) 2001-2004 European Southern Observatory
8 *
9 * This program is free software; you can redistribute it and/or modify
10 * it under the terms of the GNU General Public License as published by
11 * the Free Software Foundation; either version 2 of the License, or
12 * (at your option) any later version.
13 *
14 * This program is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 * GNU General Public License for more details.
18 *
19 * You should have received a copy of the GNU General Public License
20 * along with this program; if not, write to the Free Software
21 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
22 */
23
24 /*
25 * $Author: yjung $
26 * $Date: 2006/02/23 11:25:25 $
27 * $Revision: 1.12 $
28 * $Name: qfits-6_2_0 $
29 */
30
31 /*-----------------------------------------------------------------------------
32 Includes
33 -----------------------------------------------------------------------------*/
34
35 #include <stdio.h>
36 #include <stdlib.h>
37 #include <string.h>
38 #include <unistd.h>
39 #include <sys/types.h>
40 #include <sys/stat.h>
41
42 #include "qfits_config.h"
43 #include "qfits_image.h"
44 #include "qfits_rw.h"
45 #include "qfits_header.h"
46 #include "qfits_byteswap.h"
47 #include "qfits_tools.h"
48 #include "qfits_error.h"
49 #include "qfits_memory.h"
50 #include "qfits_std.h"
51
52 /*-----------------------------------------------------------------------------
53 Defines
54 -----------------------------------------------------------------------------*/
55
56 #define QFITSLOADERINIT_MAGIC 0xcafe
57
58 /*-----------------------------------------------------------------------------
59 Function prototypes
60 -----------------------------------------------------------------------------*/
61
62 static byte * qfits_pixdump_float(const float *, int, int);
63 static byte * qfits_pixdump_int(const int *, int, int);
64 static byte * qfits_pixdump_double(const double *, int, int);
65
66 /*----------------------------------------------------------------------------*/
67 /**
68 * @defgroup qfits_image Pixel loader for FITS images.
69 */
70 /*----------------------------------------------------------------------------*/
71 /**@{*/
72
73 /*-----------------------------------------------------------------------------
74 Function codes
75 -----------------------------------------------------------------------------*/
76
77 /*----------------------------------------------------------------------------*/
78 /**
79 @brief Dump a pixel buffer to an output FITS file in append mode.
80 @param qd qfitsdumper control object.
81 @return int 0 if Ok, -1 otherwise.
82
83 This function takes in input a qfitsdumper control object. This object
84 must be allocated beforehand and contain valid references to the data
85 to save, and how to save it.
86
87 The minimum fields to fill are:
88
89 - filename: Name of the FITS file to dump to.
90 - npix: Number of pixels in the buffer to be dumped.
91 - ptype: Type of the passed buffer (PTYPE_FLOAT, PTYPE_INT, PTYPE_DOUBLE)
92 - out_ptype: Requested FITS BITPIX for the output.
93
94 One of the following fields must point to the corresponding pixel
95 buffer:
96
97 - ibuf for an int pixel buffer (ptype=PTYPE_INT)
98 - fbuf for a float pixel buffer (ptype=PTYPE_FLOAT)
99 - dbuf for a double pixel buffer (ptype=PTYPE_DOUBLE)
100
101 This is a fairly low-level function, in the sense that it does not
102 check that the output file already contains a proper header or even
103 that the file it is appending to is indeed a FITS file. It will
104 convert the pixel buffer to the requested BITPIX type and append
105 the data to the file, without padding with zeros. See qfits_zeropad()
106 about padding.
107
108 If the given output file name is "STDOUT" (all caps), the dump
109 will be performed to stdout.
110 */
111 /*----------------------------------------------------------------------------*/
qfits_pixdump(const qfitsdumper * qd)112 int qfits_pixdump(const qfitsdumper * qd)
113 {
114 FILE * f_out;
115 byte * buf_out;
116 int buf_free;
117 int buf_sz;
118
119 /* Check inputs */
120 if (qd==NULL) return -1;
121 if (qd->filename==NULL) return -1;
122 switch (qd->ptype) {
123 case PTYPE_FLOAT:
124 if (qd->fbuf==NULL) return -1;
125 break;
126
127 case PTYPE_DOUBLE:
128 if (qd->dbuf==NULL) return -1;
129 break;
130
131 case PTYPE_INT:
132 if (qd->ibuf==NULL) return -1;
133 break;
134
135 default:
136 return -1;
137 }
138 if (qd->npix <= 0) {
139 qfits_error("Negative or NULL number of pixels specified");
140 return -1;
141 }
142
143 /*
144 * Special cases: input buffer is identical to requested format.
145 * This is only possible on big-endian machines, since FITS is
146 * big-endian only.
147 */
148 buf_out = NULL;
149 buf_free = 1;
150 #ifdef WORDS_BIGENDIAN
151 if (qd->ptype==PTYPE_FLOAT && qd->out_ptype==-32) {
152 buf_out = (byte*)qd->fbuf;
153 buf_free=0;
154 } else if (qd->ptype==PTYPE_DOUBLE && qd->out_ptype==-64) {
155 buf_out = (byte*)qd->dbuf;
156 buf_free=0;
157 } else if (qd->ptype==PTYPE_INT && qd->out_ptype==32) {
158 buf_out = (byte*)qd->ibuf;
159 buf_free=0;
160 }
161 #endif
162 buf_sz = qd->npix * BYTESPERPIXEL(qd->out_ptype);
163
164 /* General case */
165 if (buf_out==NULL) {
166 switch (qd->ptype) {
167 /* Convert buffer */
168 case PTYPE_FLOAT:
169 buf_out = qfits_pixdump_float( qd->fbuf,
170 qd->npix,
171 qd->out_ptype);
172 break;
173
174 /* Convert buffer */
175 case PTYPE_INT:
176 buf_out = qfits_pixdump_int( qd->ibuf,
177 qd->npix,
178 qd->out_ptype);
179 break;
180
181 /* Convert buffer */
182 case PTYPE_DOUBLE:
183 buf_out = qfits_pixdump_double( qd->dbuf,
184 qd->npix,
185 qd->out_ptype);
186 break;
187 }
188 }
189 if (buf_out==NULL) {
190 qfits_error("cannot dump pixel buffer");
191 return -1;
192 }
193
194 /* Dump buffer */
195 if (!strncmp(qd->filename, "STDOUT", 6)) {
196 f_out = stdout;
197 } else {
198 f_out = fopen(qd->filename, "a");
199 }
200 if (f_out==NULL) {
201 qfits_free(buf_out);
202 return -1;
203 }
204 fwrite(buf_out, buf_sz, 1, f_out);
205 if (buf_free) {
206 qfits_free(buf_out);
207 }
208 if (f_out!=stdout) {
209 fclose(f_out);
210 }
211 return 0;
212 }
213
214 /**@}*/
215
216 /*----------------------------------------------------------------------------*/
217 /**
218 @brief Convert a float pixel buffer to a byte buffer.
219 @param buf Input float buffer.
220 @param npix Number of pixels in the input buffer.
221 @param ptype Requested output BITPIX type.
222 @return 1 pointer to a newly allocated byte buffer.
223
224 This function converts the given float buffer to a buffer of bytes
225 suitable for dumping to a FITS file (i.e. big-endian, in the
226 requested pixel type). The returned pointer must be deallocated
227 using the qfits_free() function.
228 */
229 /*----------------------------------------------------------------------------*/
qfits_pixdump_float(const float * buf,int npix,int ptype)230 static byte * qfits_pixdump_float(const float * buf, int npix, int ptype)
231 {
232 byte * buf_out;
233 register byte * op;
234 int i;
235 int lpix;
236 short spix;
237 double dpix;
238
239 buf_out = qfits_malloc(npix * BYTESPERPIXEL(ptype));
240 op = buf_out;
241 switch (ptype) {
242 case 8:
243 /* Convert from float to 8 bits */
244 for (i=0; i<npix; i++) {
245 if (buf[i]>255.0) {
246 *op++ = (byte)0xff;
247 } else if (buf[i]<0.0) {
248 *op++ = (byte)0x00;
249 } else {
250 *op++ = (byte)buf[i];
251 }
252 }
253 break;
254
255 case 16:
256 /* Convert from float to 16 bits */
257 for (i=0; i<npix; i++) {
258 if (buf[i]>32767.0) {
259 *op++ = (byte)0x7f;
260 *op++ = (byte)0xff;
261 } else if (buf[i]<-32768.0) {
262 *op++ = (byte)0x80;
263 *op++ = 0x00;
264 } else {
265 spix = (short)buf[i];
266 *op++ = (spix >> 8);
267 *op++ = (spix & (byte)0xff);
268 }
269 }
270 break;
271
272 case 32:
273 /* Convert from float to 32 bits */
274 for (i=0; i<npix; i++) {
275 if (buf[i] > 2147483647.0) {
276 *op++ = (byte)0x7f;
277 *op++ = (byte)0xff;
278 *op++ = (byte)0xff;
279 *op++ = (byte)0xff;
280 } else if (buf[i]<-2147483648.0) {
281 *op++ = (byte)0x80;
282 *op++ = (byte)0x00;
283 *op++ = (byte)0x00;
284 *op++ = (byte)0x00;
285 } else {
286 lpix = (int)buf[i];
287 *op++ = (byte)(lpix >> 24);
288 *op++ = (byte)(lpix >> 16) & 0xff;
289 *op++ = (byte)(lpix >> 8 ) & 0xff;
290 *op++ = (byte)(lpix) & 0xff;
291 }
292 }
293 break;
294
295 case -32:
296 /* Convert from float to float */
297 memcpy(op, buf, npix * sizeof(float));
298 #ifndef WORDS_BIGENDIAN
299 for (i=0; i<npix; i++) {
300 qfits_swap_bytes(op, 4);
301 op++;
302 op++;
303 op++;
304 op++;
305 }
306 #endif
307 break;
308
309 case -64:
310 /* Convert from float to double */
311 for (i=0; i<npix; i++) {
312 dpix = (double)buf[i];
313 #ifndef WORDS_BIGENDIAN
314 qfits_swap_bytes(&dpix, 8);
315 #endif
316 memcpy(op, &dpix, 8);
317 op += 8;
318 }
319 break;
320
321 default:
322 qfits_error("Pixel type %i not supported yet", ptype);
323 buf_out = NULL;
324 break;
325 }
326 return buf_out;
327 }
328
329 /*----------------------------------------------------------------------------*/
330 /**
331 @brief Convert an int pixel buffer to a byte buffer.
332 @param buf Input int buffer.
333 @param npix Number of pixels in the input buffer.
334 @param ptype Requested output BITPIX type.
335 @return 1 pointer to a newly allocated byte buffer.
336
337 This function converts the given int buffer to a buffer of bytes
338 suitable for dumping to a FITS file (i.e. big-endian, in the
339 requested pixel type). The returned pointer must be deallocated
340 using the qfits_free() function.
341 */
342 /*----------------------------------------------------------------------------*/
qfits_pixdump_int(const int * buf,int npix,int ptype)343 static byte * qfits_pixdump_int(const int * buf, int npix, int ptype)
344 {
345 byte * buf_out;
346 register byte * op;
347 int i;
348
349 short spix;
350 float fpix;
351 double dpix;
352
353 buf_out = qfits_malloc(npix * BYTESPERPIXEL(ptype));
354 op = buf_out;
355 switch (ptype) {
356 case 8:
357 /* Convert from int32 to 8 bits */
358 for (i=0; i<npix; i++) {
359 if (buf[i]>255) {
360 *op++ = (byte)0xff;
361 } else if (buf[i]<0) {
362 *op++ = (byte)0x00;
363 } else {
364 *op++ = (byte)buf[i];
365 }
366 }
367 break;
368
369 case 16:
370 /* Convert from int32 to 16 bits */
371 for (i=0; i<npix; i++) {
372 if (buf[i]>32767) {
373 spix = 32767;
374 } else if (buf[i]<-32768) {
375 spix = -32768;
376 } else {
377 spix = (short)buf[i];
378 }
379 #ifndef WORDS_BIGENDIAN
380 qfits_swap_bytes(&spix, 2);
381 #endif
382 memcpy(op, &spix, 2);
383 op += 2;
384 }
385 break;
386
387 case 32:
388 /* Convert from int32 to 32 bits */
389 memcpy(op, buf, npix * sizeof(int));
390 #ifndef WORDS_BIGENDIAN
391 for (i=0; i<npix; i++) {
392 qfits_swap_bytes(op, 4);
393 op+=4;
394 }
395 #endif
396 break;
397
398 case -32:
399 /* Convert from int32 to float */
400 for (i=0; i<npix; i++) {
401 fpix = (float)buf[i];
402 #ifndef WORDS_BIGENDIAN
403 qfits_swap_bytes(&fpix, 4);
404 #endif
405 memcpy(op, &fpix, 4);
406 op += 4;
407 }
408 break;
409
410 case -64:
411 /* Convert from int32 to double */
412 for (i=0; i<npix; i++) {
413 dpix = (double)buf[i];
414 #ifndef WORDS_BIGENDIAN
415 qfits_swap_bytes(&dpix, 8);
416 #endif
417 memcpy(op, &dpix, 8);
418 op += 8;
419 }
420 break;
421
422 default:
423 qfits_error("Pixel type %i not supported yet", ptype);
424 buf_out = NULL;
425 break;
426 }
427 return buf_out;
428 }
429
430 /*----------------------------------------------------------------------------*/
431 /**
432 @brief Convert a double pixel buffer to a byte buffer.
433 @param buf Input double buffer.
434 @param npix Number of pixels in the input buffer.
435 @param ptype Requested output BITPIX type.
436 @return 1 pointer to a newly allocated byte buffer.
437
438 This function converts the given double buffer to a buffer of bytes
439 suitable for dumping to a FITS file (i.e. big-endian, in the
440 requested pixel type). The returned pointer must be deallocated
441 using the qfits_free() function.
442 */
443 /*----------------------------------------------------------------------------*/
qfits_pixdump_double(const double * buf,int npix,int ptype)444 static byte * qfits_pixdump_double(const double * buf, int npix, int ptype)
445 {
446 byte * buf_out;
447 register byte * op;
448 int i;
449
450 short spix;
451 float fpix;
452 int lpix;
453
454 buf_out = qfits_malloc(npix * BYTESPERPIXEL(ptype));
455 op = buf_out;
456 switch (ptype) {
457 case 8:
458 /* Convert from double to 8 bits */
459 for (i=0; i<npix; i++) {
460 if (buf[i]>255.0) {
461 *op++ = (byte)0xff;
462 } else if (buf[i]<0.0) {
463 *op++ = (byte)0x00;
464 } else {
465 *op++ = (byte)buf[i];
466 }
467 }
468 break;
469
470 case 16:
471 /* Convert from double to 16 bits */
472 for (i=0; i<npix; i++) {
473 if (buf[i]>32767.0) {
474 spix = 32767;
475 } else if (buf[i]<-32768.0) {
476 spix = -32768;
477 } else {
478 spix = (short)buf[i];
479 }
480 #ifndef WORDS_BIGENDIAN
481 qfits_swap_bytes(&spix, 2);
482 #endif
483 memcpy(op, &spix, 2);
484 op += 2;
485 }
486 break;
487
488 case 32:
489 /* Convert from double to 32 bits */
490 for (i=0; i<npix; i++) {
491 if (buf[i] > 2147483647.0) {
492 lpix = 2147483647;
493 } else if (buf[i] < -2147483648.0) {
494 lpix = -2147483647;
495 } else {
496 lpix = (int)buf[i];
497 }
498 #ifndef WORDS_BIGENDIAN
499 qfits_swap_bytes(&lpix, 4);
500 #endif
501 memcpy(op, &lpix, 4);
502 op += 4;
503 }
504 break;
505
506 case -32:
507 /* Convert from double to float */
508 for (i=0; i<npix; i++) {
509 fpix = (float)buf[i];
510 #ifndef WORDS_BIGENDIAN
511 qfits_swap_bytes(&fpix, 4);
512 #endif
513 memcpy(op, &fpix, 4);
514 op += 4;
515 }
516 break;
517
518 case -64:
519 /* Convert from double to double */
520 memcpy(op, buf, npix * 8);
521 #ifndef WORDS_BIGENDIAN
522 for (i=0; i<npix; i++) {
523 qfits_swap_bytes(op, 8);
524 op += 8;
525 }
526 #endif
527 break;
528
529 default:
530 qfits_error("pixel type %i not supported yet", ptype);
531 buf_out = NULL;
532 break;
533 }
534 return buf_out;
535 }
536
537 /* Test code */
538 #ifdef TESTPIXIO
qfitsloader_dump(qfitsloader * ql)539 static void qfitsloader_dump(qfitsloader * ql)
540 {
541 fprintf(stderr,
542 "file : %s\n"
543 "xtnum : %d\n"
544 "pnum : %d\n"
545 "ptype : %d\n"
546 "lx : %d\n"
547 "ly : %d\n"
548 "np : %d\n"
549 "bitpix : %d\n"
550 "seg_start : %d\n"
551 "bscale : %g\n"
552 "bzero : %g\n"
553 "ibuf : %p\n"
554 "fbuf : %p\n"
555 "dbuf : %p\n",
556 ql->filename,
557 ql->xtnum,
558 ql->pnum,
559 ql->ptype,
560 ql->lx,
561 ql->ly,
562 ql->np,
563 ql->bitpix,
564 ql->seg_start,
565 ql->bscale,
566 ql->bzero,
567 ql->ibuf,
568 ql->fbuf,
569 ql->dbuf);
570 }
571
main(int argc,char * argv[])572 int main (int argc, char * argv[])
573 {
574 qfitsloader ql;
575
576 if (argc<2) {
577 printf("use: %s <FITS>\n", argv[0]);
578 return 1;
579 }
580
581 ql.filename = argv[1];
582 ql.xtnum = 0;
583 ql.pnum = 0;
584 ql.ptype = PTYPE_FLOAT;
585
586 if (qfits_loadpix(&ql)!=0) {
587 printf("error occurred during loading: abort\n");
588 return -1;
589 }
590 qfitsloader_dump(&ql);
591 printf("pix[0]=%g\n"
592 "pix[100]=%g\n"
593 "pix[10000]=%g\n",
594 ql.fbuf[0],
595 ql.fbuf[100],
596 ql.fbuf[10000]);
597 qfits_free(ql.fbuf);
598
599 ql.ptype = PTYPE_INT;
600 if (qfits_loadpix(&ql)!=0) {
601 printf("error occurred during loading: abort\n");
602 return -1;
603 }
604 qfitsloader_dump(&ql);
605 printf("pix[0]=%d\n"
606 "pix[100]=%d\n"
607 "pix[10000]=%d\n",
608 ql.ibuf[0],
609 ql.ibuf[100],
610 ql.ibuf[10000]);
611 qfits_free(ql.ibuf);
612
613
614 ql.ptype = PTYPE_DOUBLE;
615 if (qfits_loadpix(&ql)!=0) {
616 printf("error occurred during loading: abort\n");
617 return -1;
618 }
619 qfitsloader_dump(&ql);
620 printf("pix[0]=%g\n"
621 "pix[100]=%g\n"
622 "pix[10000]=%g\n",
623 ql.dbuf[0],
624 ql.dbuf[100],
625 ql.dbuf[10000]);
626 qfits_free(ql.dbuf);
627
628 return 0;
629 }
630 #endif
631