1
2# Test routine for PDL::IO::FITS module
3
4use strict;
5
6use PDL::LiteF;
7
8use PDL::Core ':Internal'; # For howbig()
9use PDL::Config;
10
11##use PDL::Complex;  # currently not supported
12
13kill 'INT',$$  if $ENV{UNDER_DEBUGGER}; # Useful for debugging.
14
15use Test::More tests => 97;
16use Test::Exception;
17
18BEGIN {
19      use_ok( "PDL::IO::FITS" ); #1
20}
21
22require File::Spec;
23my $fs = 'File::Spec';
24sub cdir { return $fs->catdir(@_)}
25sub cfile { return $fs->catfile(@_)}
26
27my $tempd = $PDL::Config{TEMPDIR} or
28  die "TEMPDIR not found in %PDL::Config";
29my $file = cfile $tempd, "iotest$$";
30
31END {
32  unlink $file if defined $file and -e $file;
33}
34
35################ Test rfits/wfits ########################
36
37my $t = long xvals(zeroes(11,20))-5;
38
39# note: keywords are converted to uppercase
40my %hdr = ('Foo'=>'foo', 'Bar'=>42, 'NUM'=>'0123',NUMSTR=>['0123']);
41$t->sethdr(\%hdr);
42
43wfits($t, $file);
44print "#file is $file\n";
45my $t2 = rfits $file;
46
47is( sum($t->slice('0:4,:')), -sum($t2->slice('5:-1,:')),
48    "r/wfits: slice check" );				#2
49
50my $h = $t2->gethdr;
51ok( $$h{'FOO'} eq "foo" && $$h{'BAR'} == 42,
52    "header check on FOO/BAR" );			#3
53
54ok( $$h{'NUM'}+1 == 124 && $$h{'NUMSTR'} eq '0123',
55    "header check on NUM/NUMSTR" );			#4
56
57unlink $file;
58
59SKIP: {
60   eval { require Astro::FITS::Header };
61
62   skip "Astro::FITS::Header not installed", 79 if $@;
63
64########### Rudimentary table tests ################
65
66# note:
67#   the tests do not directly test the output file,
68#   instead they write out a file, read it back in, and
69#   compare to the data used to create the file.
70#   So it is more of a "self consistent" test.
71#
72sub compare_piddles ($$$) {
73    my $orig  = shift;
74    my $new   = shift;
75    my $label = shift;
76
77    TODO: {
78       local $TODO = "Need to fix alias between PDL_IND and PDL_L or PDL_LL";
79
80       is( $new->type->symbol, $orig->type->symbol, "$label has the correct type" );
81    }
82    is( $new->nelem, $orig->nelem, "  and the right number of elements" );
83    is( $new->ndims, $orig->ndims, "  and the right number of dimensions" );
84
85    my $flag;
86    if ( $orig->type() < float() ) {
87	$flag = all( $new == $orig );
88    } else {
89	$flag = all( approx( $orig, $new ) );
90    }
91    ok( $flag, "  and all the values agree" );
92}
93
94unless($PDL::Astro_FITS_Header) {
95 # Astro::FITS::Header is not present, ignore table tests
96 for(1..59){ok(1,"Test skipped (no binary table support without Astro::FITS::Header)");}
97} else { # Astro::FITS::Header exists
98
99	my $a = long( 1, 4, 9, 32 );
100	my $b = double( 2.3, 4.3, -999.0, 42 );
101	my $table = { COLA => $a, COLB => $b };
102	wfits $table, $file;
103
104	my $table2 = rfits $file;
105	unlink $file;
106
107	ok( defined $table2, "Read of table returned something" );	#5
108	is( ref($table2), "HASH", "which is a hash reference" );	#6
109	is( $$table2{tbl}, "binary", "and appears to be a binary TABLE" );#7
110
111	ok( exists $$table2{COLA} && exists $$table2{COLB}, "columns COLA and COLB exist" ); #8
112	is( $$table2{hdr}{TTYPE1}, "COLA", "column #1 is COLA" );	  #9
113	is( $$table2{hdr}{TFORM1}, "1J", "  stored as 1J" );		  #10
114	is( $$table2{hdr}{TTYPE2}, "COLB", "column #2 is COLB" );	  #11
115	is( $$table2{hdr}{TFORM2}, "1D", "  stored as 1D" );		  #12
116
117	compare_piddles $a, $$table2{COLA}, "COLA";			#13-16
118	compare_piddles $b, $$table2{COLB}, "COLB";			#17-20
119
120	$table = { BAR => $a, FOO => $b,
121		   hdr => { TTYPE1 => 'FOO', TTYPE2 => 'BAR' } };
122	$table2 = {};
123
124	wfits $table, $file;
125	$table2 = rfits $file;
126
127	ok( defined $table2 && ref($table2) eq "HASH" && $$table2{tbl} eq "binary",
128	    "Read in the second binary table" );		       #21
129	is( $$table2{hdr}{TTYPE1}, "FOO", "column #1 is FOO" );	       #22
130	is( $$table2{hdr}{TFORM1}, "1D", "  stored as 1D" );	       #23
131	is( $$table2{hdr}{TTYPE2}, "BAR", "column #2 is BAR" );	       #24
132	is( $$table2{hdr}{TFORM2}, "1J", "  stored as 1J" );	       #25
133
134	compare_piddles $a, $$table2{BAR}, "BAR";			#26-29
135	compare_piddles $b, $$table2{FOO}, "FOO";			#30-33
136
137	# try out more "exotic" data types
138
139	$a = byte(12,45,23,0);
140	$b = short(-99,100,0,32767);
141	my $c = ushort(99,32768,65535,0);
142	my $d = [ "A string", "b", "", "The last string" ];
143	my $e = float(-999.0,0,0,12.3);
144	##my $f = float(1,0,-1,2) + i * float( 0,1,2,-1 );
145	$table = {
146	       ACOL => $a, BCOL => $b, CCOL => $c, DCOL => $d, ECOL => $e,
147	##	  FCOL => $f,
148	};
149	$table2 = {};
150
151	wfits $table, $file;
152	$table2 = rfits $file;
153	#unlink $file;
154
155	ok( defined $table2 && ref($table2) eq "HASH" && $$table2{tbl} eq "binary",
156	    "Read in the third binary table" );			       #34
157	my @elem = sort keys %$table2;
158	##my @expected = sort( qw( ACOL BCOL CCOL DCOL ECOL FCOL hdr tbl ) );
159	##is ( $#elem+1, 8, "hash contains 8 elements" );
160	my @expected = sort( qw( ACOL BCOL CCOL DCOL ECOL hdr tbl ) );
161	is ( $#elem+1, 7, "hash contains 7 elements" );			#35
162	ok( eq_array( \@elem, \@expected ), "hash contains expected
163	    keys" );							#36
164
165	# convert the string array so that each element has the same length
166	# (and calculate the maximum length to use in the check below)
167	#
168	my $dlen = 0;
169	foreach my $str ( @$d ) {
170	  my $len = length($str);
171	  $dlen = $len > $dlen ? $len : $dlen;
172	}
173	foreach my $str ( @$d ) {
174	  $str .= ' ' x ($dlen-length($str));
175	}
176
177	# note that, for now, ushort data is written out as a long (Int4)
178	# instead of being written out as an Int2 using TSCALE/TZERO
179	#
180	my $i = 1;
181	foreach my $colinfo ( ( ["ACOL","1B",$a],
182				["BCOL","1I",$b],
183				["CCOL","1J",$c->long],
184				["DCOL","${dlen}A",$d],
185				["ECOL","1E",$e],
186	##			["FCOL","1M",$f]
187			      ) ) {
188	  is( $$table2{hdr}{"TTYPE$i"}, $$colinfo[0], "column $i is $$colinfo[0]" ); #37,43,49,55,58
189	  is( $$table2{hdr}{"TFORM$i"}, $$colinfo[1], "  and is stored as $$colinfo[1]" ); #38,44,50,56,59
190	  my $col = $$table2{$$colinfo[0]};
191	  if ( UNIVERSAL::isa($col,"PDL") ) {
192	    compare_piddles $col, $$colinfo[2], $$colinfo[0]; #39-42,45-48,51-54,60-63
193	  } else {
194	    # Need to somehow handle the arrays since the data read in from the
195	    # file all have 15-character length strings (or whatever the length is)
196	    #
197	    ok( eq_array($col, $$colinfo[2]),
198		"  $$colinfo[0] values agree (as an array reference)" );#57
199	  }
200	  $i++;
201	}
202}
203########### Check if r/wfits bugs are fixed ################
204
205{
206    local $| = 1;
207    my $a1 =  [1,2];
208    my $a2 = [[1,2],[1,2]];
209    my $p;
210    my $q;
211    for my $cref ( \(&byte, &short, &long, &float, &double) ) {
212        for my $a ($a1,$a2) {
213            $p = &$cref($a);
214            $p->wfits('x.fits');
215            $q = PDL->rfits('x.fits');
216	    my $flag = 1;
217            if ( ${$p->get_dataref} ne ${$q->get_dataref} ) {
218	        $flag = 0;
219	        { local $, = " ";
220		  print "\tnelem=",$p->nelem,"datatype=",$p->get_datatype,"\n";
221                  print "\tp:", unpack("c" x ($p->nelem*howbig($p->get_datatype)), ${$p->get_dataref}),"\n";
222                  print "\tq:", unpack("c" x ($q->nelem*howbig($q->get_datatype)), ${$q->get_dataref}),"\n";
223		}
224            }
225	    ok($flag,"hash reference - type check: " . &$cref ); #64-73
226        }
227    }
228    unlink 'x.fits';
229}
230
231{
232    local $| = 1;
233    my $p1= pdl  [1,2];
234    my $p2= pdl [[1,2],[1,2]];
235    my $q;
236    my @s;
237    for my $i (8,16,32,-32,-64) {
238    for my $p ($p2, $p1) {
239        $p->wfits('x.fits',$i);
240        $q = PDL->rfits('x.fits');
241        @s = $q->stats;
242	my $flag;
243	print "s=@s\n";
244        if ($s[0] == 1.5 and $s[1] < 0.7072 and $s[1]>0.577) {
245           $flag = 1;
246        } else {
247           $flag = 0;
248           print "\tBITPIX=$i, nelem=", $p->nelem, "\n";
249           print "\tbug: $s[0] == 1.5 and $s[1] == 0.5\n";
250	   { local $, = " ";
251	     print "\tp:", unpack("c8" x         $p->nelem,  ${$p->get_dataref}),"\n";
252	     print "\tq:", unpack("c" x abs($i/8*$q->nelem), ${$q->get_dataref}),"\n";
253           }
254        }
255	ok($flag,"piddle - bitpix=$i" ); #74-83
256    }
257    }
258    unlink 'x.fits';
259};
260
261}; # end of SKIP block
262
263#### Check that discontinuous data (e.g. from fftnd) get written correctly.
264#### (Sourceforge bug 3299611) it is possible to store data in a PDL non-contiguously
265#### through the C API, by manipulating dimincs; fft uses this technique, which
266#### used to hose up fits output.
267
268SKIP:{
269    eval "use PDL::FFT";
270    skip "PDL::FFT not installed", 79 if $@;
271
272    my $a = sequence(10,10,10);
273    my $ai = zeroes($a);
274    fftnd($a,$ai);
275    wfits($a,$file);
276    my $b = rfits($file);
277    ok(all($a==$b),"fftnd output (non-contiguous in memory) is written correctly");
278    unlink $file;
279}
280
281##############################
282# Check multi-HDU read/write
283
284$a = sequence(5,5);
285$b = rvals(5,5);
286
287our @aa;
288
289lives_ok { wfits([$a,$b],$file) } "wfits with multiple HDUs didn't fail";
290
291lives_ok { @aa = rfits($file) } "rfits in list context didn't fail";
292
293ok( $aa[0]->ndims == $a->ndims && all($aa[0]->shape == $a->shape), "first element has right shape");
294ok( all($aa[0] == $a), "first element reproduces written one");
295
296ok( $aa[1]->ndims == $b->ndims && all($aa[1]->shape == $b->shape), "second element has right shape");
297ok( all($aa[1] == $b), "Second element reproduces written one");
298
299unlink $file;
300
301##############################
302# Rudimentary check for longlong support
303SKIP:{
304	eval "use PDL::Types";
305	our $PDL_LL;
306    	skip "Longlong not supported",5   unless ($PDL_LL//0);
307
308	$a = rvals(longlong,7,7);
309	eval { wfits($a, $file); };
310	ok(!$@, sprintf("writing a longlong image succeeded %s",($@?"($@)":"")));
311	eval { $b = rfits($file); };
312	ok(!$@, sprintf("Reading the longlong image succeeded %s",($@?"($@)":"")));
313	ok(ref($b->hdr) eq "HASH", "Reading the longlong image produced a PDL with a hash header");
314	ok($b->hdr->{BITPIX} == 64, "BITPIX value was correct");
315	ok(all($b==$a),"The new image matches the old one (longlong)");
316	unlink $file;
317}
318
319###############################
320# Check that tilde expansion works
321my $tildefile = cfile('~',"PDL-IO-FITS-test_$$.fits");
322lives_ok { sequence(3,5,2)->wfits($tildefile) } "wfits tilde expansion didn't fail";
323lives_ok { rfits($tildefile) } "rfits tilde expansion didn't fail";
324$tildefile =~ s/^(~)/glob($1)/e; #use the same trick as in FITS.pm to resolve this filename.
325unlink($tildefile) or warn "Could not delete $tildefile: $!\n"; #clean up.
326
3271;
328