I am trying to translate a python code into a C code. So I was trying to translate that line:
import healpy as hp
OldTobySensitivMap = hp.read_map('file.fits', dtype=numpy.float64)
This opens a FITS file using the HEALPix projection. So in C I am using the HEALPix library:
#include "chealpix.h"
#include "fitsio.h"
#include <string.h>
int main(int argc, char **argv){
char coordsys[]="G";
char order[]="NESTED";
long nside=0;
float *map;
map=read_healpix_map("file.fits",&nside,coordsys,order);
printf("nside:%ld, coordsys:%s, order:%s",nside,coordsys,order);
return(0);
}
However, this code returns an error :
WARNING: Could not find COORDSYS keyword in in file file.fits
*** stack smashing detected ***: terminated
Program received signal SIGABRT, Aborted.
GDB debugger output is:
Program received signal SIGABRT, Aborted.
__GI_raise (sig=sig@entry=6) at ../sysdeps/unix/sysv/linux/raise.c:50
50 ../sysdeps/unix/sysv/linux/raise.c: No such file or directory.
(gdb) bt
#0 __GI_raise (sig=sig@entry=6) at ../sysdeps/unix/sysv/linux/raise.c:50
#1 0x00007ffff76f5859 in __GI_abort () at abort.c:79
#2 0x00007ffff77603ee in __libc_message (action=action@entry=do_abort, fmt=fmt@entry=0x7ffff788a07c "*** %s ***: terminated\n") at ../sysdeps/posix/libc_fatal.c:155
#3 0x00007ffff7802b4a in __GI___fortify_fail (msg=msg@entry=0x7ffff788a064 "stack smashing detected") at fortify_fail.c:26
#4 0x00007ffff7802b16 in __stack_chk_fail () at stack_chk_fail.c:24
#5 0x00005555555556af in main (argc=1, argv=0x7fffffffdd18) at young_pulsar_pop.c:60
Increasing the size of the strings such that:
char coordsys[2024]="G";
char order[2024]="NESTED";
Works, but then the output is nside:512, coordsys:, order:RING, so I do not really understand what is causing the stack error.
Like it was said by an other user, that HEALPix interface seem to be quite unsafe for C, and the documentation is really poor. I am therefore wondering if I should not simply embed the python code in my C programm.
If you check the related documentation :
read_healpix_map
This routine reads a full sky HEALPix map from a FITS file
Location in HEALPix directory tree: src/C/subs/read_healpix_map.c
FORMAT float *read_healpix_map(char *infile, long *nside, char *coordsys, char *ordering)
ARGUMENTS
| name&dimensionality | kind | in/out | description | read_healpix_map float OUT array containing the map read from the file infile char IN FITS file containing a full sky to be read nside long OUT HEALPix resolution parameter of the map coordsys char OUT astronomical coordinate system of pixelisation (either 'C', 'E' or 'G' standing respectively for Celestial=equatorial, Ecliptic or Galactic) ordering char OUT HEALPix pixel ordering (either 'RING' or 'NESTED')
What we see is that coordsys and ordering are output parameters for this function. That means this function expect char buffers with a defined size to store the related output string. As described in this documentation, the provided buffer for coordsys should be, for ASCII string a char [2] type (storing the coordinate system in one char plus the NULL byte) and ordering should be a char [7] type (the longest string is 'NESTED', 6 char plus the NULL byte).
With char coordsys[]="G"; and char order[]="NESTED";, you're defining the string in the data section of your program. As you use the longest string for order, you should then have sufficient space for the both strings.
The read_healpix_map documentation isn't accurate. The function will store a TSTRING object in its coordsys and ordering parameters, as its source calls the fits_read_key function which store a TSRING object in theses parameters.
So that's the reason why using 'big size' buffers works.
The best you can do to avoid any troubles is to directly initialize your buffers coordsys and order with a TSTRING object for "G" and "NESTED". In the fitsio.h header documentation, the max length for a keyword is defined with FLEN_VALUE. So i suggest you to use:
char coordsys[FLEN_VALUE]="G";
char order[FLEN_VALUE]="NESTED";