1# -*-perl-*- 2# 3# test bad value handling in PDL 4# - as it's a compile-time option we 5# skip unless $PDL::Config{WITH_BADVAL} 6# 7 8use strict; 9use Test::More; 10 11# although approx() caches the tolerance value, we 12# use it in every call just to document things 13# 14use constant ABSTOL => 1.0e-4; 15 16use File::Temp qw( tempfile ); 17my $fname; 18{ 19 local $^W = 0; 20 (undef, $fname) = tempfile( 'delmeXXXXX', SUFFIX => '.fits', OPEN => 0 ); 21} 22 23END { 24 unlink $fname if -e $fname; 25} 26 27use PDL::LiteF; 28$| = 1; 29 30use PDL::Config; 31if ( $PDL::Config{WITH_BADVAL} ) { 32 plan tests => 82; 33} else { 34 # reduced testing 35 plan tests => 10; 36 37 my $a = pdl(1,2,3); 38 is( $a->badflag(), 0 ); # 1 39 40 my $b = pdl(4,5,6); 41 my $c = $a + $b; 42 is( $c->badflag(), 0 ); # 2 43 is( $c->sum(), 21 ); # 3 44 45 # can not set the bad flag 46 $a->badflag(1); 47 is( $a->badflag(), 0 ); 48 49 # and piddles do not have a bad value 50 ok( ! defined $a->badvalue ); 51 52 # can not change a piddle to include bad values 53 ok( all( ($a->setbadif( $a == 2 ) - pdl(1,2,3)) == 0 ) ); 54 55 $a = ones(3,2,4); 56 $b = zeroes(2,4); 57 $c = ones(2,4) * 3; 58 is( $a->nbad, 0 ); 59 is( $a->ngood, 24 ); 60 ok( all( ($a->nbadover - $b) == 0 ) ); 61 ok( all( ($a->ngoodover - $c) == 0 ) ); 62 63 exit; 64} 65 66my $usenan = $PDL::Config{BADVAL_USENAN} || 0; 67my $perpdl = $PDL::Config{BADVAL_PER_PDL} || 0; 68 69# check default behaviour (ie no bad data) 70# - probably overkill 71# 72my $a = pdl(1,2,3); 73is( $a->badflag(), 0, "no badflag" ); 74 75my $b = pdl(4,5,6); 76my $c = $a + $b; 77is( $c->badflag(), 0, "badflag not set in a copy" ); 78is( $c->sum(), 21, "sum() works on non bad-flag piddles" ); 79 80# is the flag propagated? 81$a->badflag(1); 82ok( $a->badflag(), "bad flag is now set" ); 83 84$c = $a + $b; 85ok( $c->badflag(), "bad flag is propagated" ); 86is( $c->sum(), 21, "sum is still 21 with badflag set" ); 87 88$a->badflag(0); 89$b->badflag(1); 90$c = $a + $b; 91ok( $c->badflag(), "badflag propagates on rhs of 'a+b'" ); 92 93# how about copies/vaffines/whatever 94$a = rvals( long, 7, 7, {Centre=>[2,2]} ); 95$b = $a; 96is( $b->badflag, 0, "badflag not set in a copy" ); 97 98$a->badflag(1); 99$b = $a; 100ok( $b->badflag, "badflag is now set in a copy" ); 101 102$a->badflag(0); 103$b = $a->slice('2:5,3:4'); 104$c = $b->slice('0:1,(0)'); 105is( $b->badflag, 0, "slice handling okay with no badflag" ); 106 107$a->badflag(1); 108 109my $i = "Type: %T Dim: %-15D State: %5S Dataflow: %F"; 110print "Info: a = ", $a->info($i), "\n"; 111print "Info: b = ", $b->info($i), "\n"; 112print "Info: c = ", $b->info($i), "\n"; 113 114# let's check that it gets through to a child of a child 115ok( $c->badflag, "badflag propagated throufh to a child" ); 116 117# can we change bad values 118is( byte->badvalue, byte->orig_badvalue, "byte bad value is set to the default value" ); 119byte->badvalue(23); 120is( byte->badvalue, 23, "changed bad value for byte" ); 121byte->badvalue( byte->orig_badvalue ); 122 123# check setbadat() 124$a = pdl(1,2,3,4,5); 125$a->setbadat(2); 126is( PDL::Core::string($a), "[1 2 BAD 4 5]", "setbadat worked" ); 127 128# now check that badvalue() changes the piddle 129# (only for integer types) 130$a = convert($a,ushort); 131my $badval = $a->badvalue; 132$a->badvalue(44); 133is( PDL::Core::string($a), "[1 2 BAD 4 5]", "changed badvalue" ); 134$a->badflag(0); 135is( PDL::Core::string($a), "[1 2 44 4 5]", "can remove the bad value setting" ); 136 137# restore the bad value 138$a->badvalue($badval); 139 140$a = byte(1,2,3); 141$b = byte(1,byte->badvalue,3); 142$a->badflag(1); 143$b->badflag(1); 144 145# does string work? 146# (this has implicitly been tested just above) 147# 148is( PDL::Core::string($b), "[1 BAD 3]", "can convert bad values to a string" ); 149 150# does addition work 151$c = $a + $b; 152is( sum($c), 8, "addition propagates the bad value" ); 153 154# does conversion of bad types work 155$c = float($b); 156ok( $c->badflag, "type conversion retains bad flag" ); 157is( PDL::Core::string($c), "[1 BAD 3]", " and the value" ); 158is( sum($c), 4, " and the sum" ); 159 160$a = byte(1,2,byte->badvalue,byte->badvalue,5,6,byte->badvalue,8,9); 161$a->badflag(1); 162 163is( PDL::Core::string($a->isbad), "[0 0 1 1 0 0 1 0 0]", "isbad() works" ); 164is( PDL::Core::string($a->isgood), "[1 1 0 0 1 1 0 1 1]", "isgood() works" ); 165 166is( $a->nbad, 3, "nbad() works" ); 167is( $a->ngood, 6, "ngood() works" ); 168 169$a = byte( [255,255], [0,255], [0,0] ); 170$a->badflag(1); 171 172is( PDL::Core::string($a->nbadover), "[2 1 0]", "nbadover() works" ); 173is( PDL::Core::string($a->ngoodover), "[0 1 2]", "ngoodover() works" ); 174 175# check dataflow (or vaffine or whatever it's called) 176$a = byte( [1,2,byte->badvalue,4,5], [byte->badvalue,0,1,2,byte->badvalue] ); 177$a->badflag(1); 178$b = $a->slice(',(1)'); 179is( sum($b), 3, "sum of slice works" ); 180$b++; 181is( PDL::Core::string($a), 182 "\n[\n [ 1 2 BAD 4 5]\n [BAD 1 2 3 BAD]\n]\n", 183 "inplace addition of slice flows back to parent" 184 ); 185 186$a = byte->badvalue * ones(byte,3,2); 187is( $a->get_datatype, 0, "datatype remains a byte" ); 188$a->badflag(1); 189is( PDL::Core::string( PDL::zcover($a) ), "[BAD BAD]", "zcover() okay" ); 190$a->set(1,1,1); 191$a->set(2,1,1); 192is( PDL::Core::string( PDL::zcover($a) ), "[BAD 0]", " and still okay" ); 193 194# 255 is the default bad value for a byte array 195# 196$a = byte(1,2,255,4,5); 197is( $a->median, 4, "median() works on good piddle" ); 198$a->badflag(1); 199is( $a->median, 3, "median() works on bad biddle" ); 200 201# as random() creates numbers between 0 and 1 it won't 202# accidentally create a bad value by chance (the default 203# bad value for a double is either a very negative 204# number or NaN). 205# 206$a = random(20); 207$a->badflag(1); 208is( $a->check_badflag, 0, "check_badflag did not find a bad value" ); 209 210$i = "Type: %T Dim: %-15D State: %5S Dataflow: %F"; 211 212# check out stats, since it uses several routines 213# and setbadif 214$a = pdl( qw(42 47 98 13 22 96 74 41 79 76 96 3 32 76 25 59 5 96 32 6) ); 215$b = $a->setbadif( $a < 20 ); 216my @s = $b->stats(); 217ok( approx( $s[0], 61.9375, ABSTOL ), "setbadif/stats test 1" ); 218ok( approx( $s[1], 27.6079, ABSTOL ), "setbadif/stats test 2" ); 219is( $s[2], 66.5, "setbadif/stats test 3" ); 220is( $s[3], 22, "setbadif/stats test 4" ); 221is( $s[4], 98, "setbadif/stats test 5" ); 222ok( approx( $s[6], 26.7312, ABSTOL ), "setbadif/stats test 6" ); 223 224# how about setbadtoval (was replacebad) 225$a = $b->setbadtoval(20) - pdl(qw(42 47 98 20 22 96 74 41 79 76 96 20 32 76 25 59 20 96 32 20)); 226ok( all($a == 0), "setbadtoval() worked" ); 227 228# and inplace? 229$a = pdl( qw(42 47 98 13 22 96 74 41 79 76 96 3 32 76 25 59 5 96 32 6) ); 230$b = $a->setbadif( $a < 20 ); 231$b->inplace->setbadtoval(20); 232$a = $b - pdl(qw(42 47 98 20 22 96 74 41 79 76 96 20 32 76 25 59 20 96 32 20)); 233ok( all($a == 0), " and inplace" ); 234 235# ditto for copybad 236$a = pdl( qw(42 47 98 13 22 96 74 41 79 76 96 3 32 76 25 59 5 96 32 6) ); 237$b = $a->setbadif( $a < 20 ); 238$c = copybad( $a, $b ); 239is( PDL::Core::string( $c->isbad ), 240 "[0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 1]", 241 "isbad() worked" ); 242 243$a = pdl( qw(42 47 98 13 22 96 74 41 79 76 96 3 32 76 25 59 5 96 32 6) ); 244$b = $a->setbadif( $a < 20 ); 245$a->inplace->copybad( $b ); 246is( PDL::Core::string( $a->isbad ), 247 "[0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 1]", 248 " and inplace" ); 249 250## $a->inplace->setbadif( $a % 2 ) does NOT work because 251## ($a % 2) is performed inplace - ie the flag is set for 252## that function 253# 254##$a = sequence(3,3); 255##$a->inplace->setbadif( $a % 2 ); 256###$a = $a->setbadif( $a % 2 ); # for when not bothered about inplace 257##ok( PDL::Core::string( $a->clump(-1) ), 258## "[0 BAD 2 BAD 4 BAD 6 BAD 8]" ); # 259 260## look at propagation of bad flag using inplace routines... 261$a = sequence( byte, 2, 3 ); 262$a = $a->setbadif( $a == 3 ); 263$b = $a->slice("(1),:"); 264$a->inplace->setbadtoval(3); 265is( $b->badflag, 0, "badflag cleared using inplace setbadtoval()" ); 266 267$a = sequence( byte, 2, 3 ); 268$b = $a->slice("(1),:"); 269my $mask = sequence( byte, 2, 3 ); 270$mask = $mask->setbadif( ($mask % 3) == 2 ); 271print "a,b == ", $a->badflag, ",", $b->badflag, "\n"; 272$a->inplace->copybad( $mask ); 273print "a,b == ", $a->badflag, ",", $b->badflag, "\n"; 274print "$a $b\n"; 275is( $b->badflag, 1, "badflag propagated using inplace copybad()" ); 276 277# test some of the qsort functions 278$a = pdl( qw(42 47 98 13 22 96 74 41 79 76 96 3 32 76 25 59 5 96 32 6) ); 279$b = $a->setbadif( $a < 20 ); 280my $ix = qsorti( $b ); 281is( PDL::Core::string( $b->index($ix) ), 282 "[22 25 32 32 41 42 47 59 74 76 76 79 96 96 96 98 BAD BAD BAD BAD]", 283 "qsorti() okay" 284 ); # 285 286# check comparison/bit operators in ops.pd 287 288$a = pdl( 2, 4, double->badvalue ); 289$a->badflag(1); 290$b = abs( $a - pdl(2.001,3.9999,234e23) ) > 0.01; 291is( PDL::Core::string( $b ), "[0 0 BAD]", "abs() and >" ); 292 293$b = byte(1,2,byte->badvalue,4); 294$b->badflag(1); 295is( PDL::Core::string( $b << 2 ), "[4 8 BAD 16]", "<<" ); 296 297$a = pdl([1,2,3]); 298$a->badflag(1); 299$b = $a->assgn; 300is( $b->badflag, 1, "assgn propagated badflag"); 301$a->badflag(0); 302is( $b->badflag, 1, "assgn is not a deep copy for the badflag"); 303 304# check that at and sclr return the correct values 305TODO: { 306 $a = pdl q[BAD]; 307 local $TODO = 'check that at and sclr return correct values and the same'; 308 309 is( PDL::Core::string($a), 'BAD', 'can convert PDL to string' ); 310 is( $a->at, 'BAD', 'at() returns BAD for a bad value' ); 311 is( $a->sclr, 'BAD', 'sclr() returns BAD for a bad value' ); 312} 313 314# quick look at math.pd 315use PDL::Math; 316 317$a = pdl(0.5,double->badvalue,0); 318$a->badflag(1); 319$b = bessj0($a); 320is( PDL::Core::string( isbad($b) ), "[0 1 0]", "bessj0()" ); 321 322$a = pdl(double->badvalue,0.8); 323$a->badflag(1); 324$b = bessjn($a,3); # thread over n() 325is( PDL::Core::string( isbad($b) ), "[1 0]", "thread over bessjn()" ); 326ok( abs($b->at(1)-0.010) < 0.001 ); 327 328$a = pdl( 0.01, 0.0 ); 329$a->badflag(1); 330ok( all( abs(erfi($a)-pdl(0.00886,0)) < 0.001 ), "erfi()" ); 331 332# I haven't changed rotate, but it should work anyway 333$a = byte( 0, 1, 2, 4, 5 ); 334$a->setbadat(2); 335is( PDL::Core::string( $a->rotate(2) ), "[4 5 0 1 BAD]", "rotate()" ); 336 337# check norm 338$a = float( 2, 0, 2, 2 )->setvaltobad(0.0); 339$b = $a->norm; 340$c = $a/sqrt(sum($a*$a)); 341ok( all( approx( $b, $c, ABSTOL ) ), "norm()" ); 342 343# Image2D 344my $ans = pdl( 345 [ 3, 7, 11, 21, 27, 33, 39, 45, 51, 27], 346 [ 3, 5, 13, 21, 27, 33, 39, 45, 51, 27], 347 [ 3, 9, 15, 21, 27, 33, 39, 45, 51, 27] 348); 349 350$a = xvals zeroes 10,3; 351$a->setbadat(2,1); 352 353$b = pdl [1,2],[2,1]; 354 355use PDL::Image2D; 356$c = conv2d($a, $b); 357 358is( int(at(sum($c-$ans))), 0, "conv2d()" ); 359 360$a = zeroes(5,5); 361$a->badflag(1); 362my $t = $a->slice("1:3,1:3"); 363$t .= ones(3,3); 364$a->setbadat(2,2); 365 366$b = sequence(3,3); 367$ans = pdl ( [0,0,0,0,0],[0,0,2,0,0],[0,1,5,2,0],[0,0,4,0,0],[0,0,0,0,0]); 368is( int(at(sum(med2d($a,$b)-$ans))), 0, "med2d()" ); 369 370# propagation of badflag using inplace ops (ops.pd) 371 372# test biop fns 373$a = sequence(3,3); 374$c = $a->slice(',(1)'); 375$b = $a->setbadif( $a % 2 ); 376$a->inplace->plus($b,0); 377print $a; 378print "$c\n"; 379is( PDL::Core::string($c), "[BAD 8 BAD]", "inplace biop - plus()" ); 380 381# test bifunc fns 382$a = sequence(3,3); 383$c = $a->slice(',(1)'); 384$b = $a->setbadif( $a % 3 != 0 ); 385$a->inplace->power($b,0); 386print $a; 387print "$c\n"; 388is( PDL::Core::string($c), "[27 BAD BAD]", "inplace bifunc - power()" ); 389 390# test histogram (using hist) 391$a = pdl( qw/1 2 3 4 5 4 3 2 2 1/ ); 392$a->setbadat(1); 393$b = hist $a, 0, 6, 1; 394print "values: $a\n"; 395print "histogram: $b\n"; 396is( PDL::Core::string($b), "[0 2 2 2 2 1]", "hist()" ); 397 398#$b = $a->isfinite; 399#print "isfinite(A): datatype = [",$b->get_datatype,"]\n"; 400 401$a->inplace->isfinite; 402#print "A: datatype = [",$a->get_datatype,"]\n"; 403is( PDL::Core::string($a), "[1 0 1 1 1 1 1 1 1 1]", "isfinite()" ); 404#print "A: datatype = [",$a->get_datatype,"]\n"; 405 406# histogram2d 407$a = long(1,1,1,2,2); 408$b = long(2,1,1,1,1); 409$b->setbadat(0); 410my @c = ( 1,0,3 ); 411$c = histogram2d($a,$b,@c,@c); 412is( PDL::Core::string($c->clump(-1)), 413 "[0 0 0 0 2 2 0 0 0]", 414 "histogram2d()" ); 415 416# weird propagation of bad values 417# - or is it? 418# 419#$a = sequence( byte, 2, 3 ); 420#$a = $a->setbadif( $a == 3 ); 421#$b = $a->slice("(1),:"); 422#$a .= $a->setbadtoval(3); 423#ok( $a->badflag, 0 ); # this fails 424#ok( $b->badflag, 0 ); # as does this 425 426# badmask: inplace 427$a = sequence(5); 428$a->setbadat(2); 429$a->inplace->badmask(0); 430is( PDL::Core::string($a), "[0 1 0 3 4]", "inplace badmask()" ); 431 432# setvaltobad 433$a = sequence(10) % 4; 434$a->inplace->setvaltobad( 1 ); 435like( PDL::Core::string( $a->clump(-1) ), 436 qr{^\[-?0 BAD 2 3 -?0 BAD 2 3 -?0 BAD]$}, "inplace setvaltobad()" ); 437 438# check setvaltobad for non-double piddles 439my $fa = pdl( float, 1..4) / 3; 440my $da = pdl( double, 1..4) / 3; 441ok( all($fa->setvaltobad(2/3)->isbad == $da->setvaltobad(2/3)->isbad), "setvaltobad for float piddle"); 442 443# simple test for setnantobad 444# - could have a 1D FITS image containing 445# NaN's and then a simple version of rfits 446# (can't use rfits as does conversion!) 447$a->inplace->setnantobad; 448like( PDL::Core::string( $a->clump(-1) ), 449 qr{^\[-?0 BAD 2 3 -?0 BAD 2 3 -?0 BAD]$}, "inplace setnantobad()" ); 450 451# test r/wfits 452use PDL::IO::FITS; 453 454$a = sequence(10)->setbadat(0); 455print "Writing to fits: $a type = (", $a->get_datatype, ")\n"; 456$a->wfits($fname); 457$b = rfits($fname); 458print "Read from fits: $b type = (", $b->get_datatype, ")\n"; 459 460ok( $b->slice('0:0')->isbad, "rfits/wfits propagated bad flag" ); 461ok( sum(abs($a-$b)) < 1.0e-5, " and values" ); 462 463# now force to integer 464$a->wfits($fname,16); 465$b = rfits($fname); 466print "BITPIX 16: datatype == ", $b->get_datatype, " badvalue == ", $b->badvalue(), "\n"; 467ok( $b->slice('0:0')->isbad, "wfits coerced bad flag with integer datatype" ); 468ok( sum(abs(convert($a,short)-$b)) < 1.0e-5, " and the values" ); 469 470# check that we can change the value used to represent 471# missing elements for floating points (earlier tests only did integer types) 472# IF we are not using NaN's 473# 474SKIP: { 475 skip( "Skipped: test only valid when not using NaN's as bad values", 2 ) 476 if $usenan; 477 478 # perhaps should check that the value can't be changed when NaN's are 479 # being used. 480 # 481 482 is( float->badvalue, float->orig_badvalue, "default bad value for floats matches" ); 483 is( float->badvalue(23), 23, "changed floating-point bad value" ); 484 float->badvalue( float->orig_badvalue ); 485} 486 487SKIP: { 488 489 skip ("Skipped: test only valid when enabling bad values per pdl", 3) 490 unless $perpdl; 491 492 $a = sequence(4); 493 $a->badvalue(3); 494 $a->badflag(1); 495 $b = $a->slice('2:3'); 496 is( $b->badvalue, 3, "can propagate per-piddle bad value"); 497 is( $b->sum, 2, "and the propagated value is recognised as bad"); 498 499 $a = sequence(4); 500 is ($a->badvalue, double->orig_badvalue, "no long-term affects of per-piddle changes [1]"); 501 502} 503 504SKIP: { 505 skip ("Skipped: test not valid if per-piddle bad values are used", 1) 506 if $perpdl; 507 508 $a = double(4); 509 double->badvalue(3); 510 is($a->badvalue, double->badvalue, "no long-term affects of per-piddle changes [2]"); 511 double->badvalue(double->orig_badvalue); 512 513} 514 515# At the moment we do not allow per-piddle bad values 516# and the use of NaN's. 517#TODO: { 518# local $TODO = "Need to work out whan NaN and per-piddle bad values means"; 519# is (0, 1); 520#} 521 522# end 523