diff --git a/Libtmp/Slatec/t/slatec.t b/Libtmp/Slatec/t/slatec.t
index 82c53b4c6..d1c0fc335 100644
--- a/Libtmp/Slatec/t/slatec.t
+++ b/Libtmp/Slatec/t/slatec.t
@@ -2,36 +2,27 @@ use strict;
 use warnings;
 use PDL::LiteF;
 use Test::More;
+use Test::PDL;
 use PDL::Slatec;
 use PDL::MatrixOps qw(identity);
 
 kill 'INT', $$ if $ENV{UNDER_DEBUGGER}; # Useful for debugging.
 
-sub tapprox {
-	my($x,$y,$c,$d) = @_;
-	$c = abs($x-$y);
-	$d = max($c);
-#	print "APR: $x,$y,$c,$d;\n";
-	$d < 0.001;
-}
-
 my $mat = pdl [1,0.1],[0.1,2];
 
 my ($eigvals,$eigvecs) = eigsys($mat);
 
 ## print STDERR $eigvecs,$eigvals,"\n";
 
-ok(tapprox($eigvals,pdl(0.9901,2.009)));
-ok(!tapprox($eigvals,pdl(0.99,2.5)));
-
-ok(tapprox($eigvecs,pdl([0.995,-0.0985],[0.0985,0.995])));
+is_pdl $eigvals,float(0.9901,2.009), {atol=>1e-3};
+is_pdl $eigvecs,float([0.995,-0.0985],[0.0985,0.995]), {atol=>1e-3};
 
 $mat = pdl [2,3],[4,5];
 
 my $inv = matinv($mat);
 
 my $uni=scalar $mat x $inv;
-ok(tapprox($uni,identity(2)));
+is_pdl $uni,identity(2);
 
 eval {matinv(identity(2)->dummy(-1,2))};
 is $@, '', 'matinv can broadcast';
@@ -39,8 +30,8 @@ is $@, '', 'matinv can broadcast';
 my $det = $mat->det;
 my $deti = $inv->det;
 
-ok(tapprox($det,-2));
-ok(tapprox($deti,-0.5));
+is_pdl $det, pdl(-2);
+is_pdl $deti, pdl(-0.5);
 
 # Now do the polynomial fitting tests
 
@@ -133,38 +124,34 @@ $xx    = pdl(12,4,6.25,1.5); # Ask for multiple positions at once
 
 ($yfit, $yp) = polyvalue($ndeg, $nder, $xx, $a1);
 
-## print STDERR "At $xx is $yfit and $yp\n";
-
 # Simple test of expected value                                               
 ok(int($yfit->at(1)) == 15);            
 
 # test the PCHIP stuff
 
 ## Test: chim/chic
-#
 $x = float( 3 .. 10 );
 my $f = $x*$x*$x + 425.42352;
-my $answer = 3.0*$x*$x; 
+my $answer = 3*$x*$x;
 
 my ( $d, $err ) = chim( float($x), float($f) );
-
-ok(($err->getndims==0) & ($err->sum == 0));
+is_pdl $err, longlong 0;
 
 # don't check the first and last elements, as expect the
 # error to be largest there
 # value of 5% comes from tests on linux and solaris machines
-ok(all( slice( abs(($d - $answer)/$answer), '1:-2' ) < 0.05 ) );
+is_pdl +(map $_->slice('1:-2'), $d, $answer), {rtol=>0.05};
 
 # compare the results of chic
 my $d2 = $f->zeroes;
 chic( pdl([0, 0]), pdl([0, 0]), 1, $x, $f, $d2, my $err2=null );
-ok(($err2->getndims==0) & ($err2->sum == 0));
-ok(all( abs($d2 - $d) < 0.02 ) );
+is_pdl $err2, longlong 0;
+is_pdl $d2, $d->double, {atol=>2e-2};
 
 ## Test: chsp
 chsp( pdl([0, 0]), pdl([0, 0]), $x, $f, my $d3=null, my $err3=null );
-ok(($err3->getndims==0) & ($err3->sum == 0));
-ok(all( abs($d3 - $d) < 2 ) );
+is_pdl $err3, longlong 0;
+is_pdl $d3, $d->double, {atol=>2};
 
 ## Test: chfd/chfe
 #
@@ -175,33 +162,27 @@ my ( $fe, $de );
 ok(($err->getndims==0) & ($err->sum == 0));
 
 $answer = $xe*$xe*$xe + 425.42352;
-ok(all( abs(($fe - $answer)/$answer) < 1.0e-5 ) );
+is_pdl $fe, $answer, {rtol=>1e-5};
 
 $answer = 3.0*$xe*$xe;
-ok(all( abs(($de - $answer)/$answer) < 0.02 ) );
+is_pdl $de, $answer, {rtol=>2e-2};
 
 ( $fe, $err ) = chfe( $x, $f, $d, 1, $xe );
-
-ok(($err->getndims==0) & ($err->sum == 0));
-
-$answer = $xe*$xe*$xe + 425.42352;
-ok(all( abs(($fe - $answer)/$answer) < 1.0e-5 ) );
+is_pdl $fe, $xe*$xe*$xe + 425.42352, {rtol=>1e-3};
+is_pdl $err, longlong 0;
 
 # Test: chcm
 #
 $x   = float( 1, 2, 3, 5, 6, 7 );
 $f   = float( 1, 2, 3, 4, 3, 4 );
-my $ans = long(  1, 1, 1, -1, 1, 2 );
+my $ans = longlong(  1, 1, 1, -1, 1, 2 );
 
 ( $d, $err ) = chim($x, $f);
-ok(($err->getndims==0) & ($err->sum == 2)); # 2 switches in monotonicity
-
-my $ismon;
-( $ismon, $err ) = chcm($x, $f, $d, 1);
+is_pdl $err, longlong 2;
 
-ok(($err->getndims==0) & ($err->sum == 0));
-ok($ismon->get_datatype == longlong->enum) or diag $ismon->type;
-ok(tapprox($ismon,$ans)) or diag $ismon, "\nexpected:\n", $ans;
+( my $ismon, $err ) = chcm($x, $f, $d, 1);
+is_pdl $err, longlong 0;
+is_pdl $ismon,$ans;
 
 ## Test: chia
 #
@@ -211,17 +192,15 @@ $f = $x * $x;
 
 $ans = pdl( 9.0**3, (8.0**3-1.0**3) ) / 3.0;
 ( my $int, $err ) = chia($x, $f, $d, my $skip=zeroes(2), pdl(0.0,1.0), pdl(9.0,8.0));
-ok(all($err == 0));
-ok(all( abs($int-$ans) < 0.04 ) );
+is_pdl $err, longlong 0,0;
+is_pdl $int, $ans, {atol=>4e-2};
 
 my $hi = pdl( $x->at(9), $x->at(7) );
 my $lo = pdl( $x->at(0), $x->at(1) );
 $ans = ($hi**3 - $lo**3) / 3;
 ( $int, $err ) = chid( $x, $f, $d, $skip=zeroes(2), pdl(0,1), pdl(9,7) );
-ok(all($err == 0));
-ok(all( abs($int-$ans) < 0.06 ) );
-## print STDERR "int=$int; ans=$ans; int-ans=".($int-$ans)."\n";
-## print STDERR "ref ans=".(ref $ans)."\n";
+is_pdl $err, longlong 0,0;
+is_pdl $int, $ans, {atol=>6e-2};
 
 ## Test: chbs - note, only tests that it runs successfully
 my $nknots = 0;
@@ -231,7 +210,7 @@ my $ndim = PDL->null;
 my $kord = PDL->null;
 $err = PDL->null;
 chbs( $x, $f, $d, 0, $nknots, $t, $bcoef, $ndim, $kord, $err );
-ok(all($err == 0));
+is_pdl $err, longlong 0;
 
 ## Test: bvalu - note, only tests that it runs successfully
 my $x_slice = $x->slice('0:-2'); # because calling with last value is out of range
@@ -244,21 +223,21 @@ gefa(my $lu=$A->copy, my $ipiv=null, my $info=null);
 gesl($lu, $ipiv, $x=$B->transpose->copy, 1); # 1 = do transpose because Fortran
 $x = $x->inplace->transpose;
 my $got = $A x $x;
-ok tapprox $got, $B or diag "got: $got";
+is_pdl $got, $B;
 
 {
 my $pa = pdl(float,1,-1,1,-1); # even number
 my ($az, $x, $y) = PDL::Slatec::fft($pa);
-ok all approx $az, 0;
-ok all approx $x, pdl "[0 1 0 0]";
-ok all approx $y, pdl "[0 0 0 0]";
-ok all approx PDL::Slatec::rfft($az, $x, $y), $pa;
+is_pdl $az, float 0;
+is_pdl $x, float "[0 1 0 0]";
+is_pdl $y, float "[0 0 0 0]";
+is_pdl PDL::Slatec::rfft($az, $x, $y), $pa;
 $pa = pdl(float,1,-1,1,-1,1); # odd number
 ($az, $x, $y) = PDL::Slatec::fft($pa);
-ok all approx $az, 0.2;
-ok all approx $x, pdl "[0.4 0.4 0 0 0]";
-ok all approx $y, pdl "[-0.2906170 -1.231073 0 0 0]";
-ok all approx PDL::Slatec::rfft($az, $x, $y), $pa;
+is_pdl $az, float 0.2;
+is_pdl $x, float "[0.4 0.4 0 0 0]";
+is_pdl $y, float "[-0.2906170 -1.231073 0 0 0]";
+is_pdl PDL::Slatec::rfft($az, $x, $y), $pa;
 }
 
 done_testing;
diff --git a/Libtmp/Transform-Proj4/t/gis_proj.t b/Libtmp/Transform-Proj4/t/gis_proj.t
index 1960f1386..8959156b9 100644
--- a/Libtmp/Transform-Proj4/t/gis_proj.t
+++ b/Libtmp/Transform-Proj4/t/gis_proj.t
@@ -2,6 +2,7 @@ use strict;
 use warnings;
 use PDL::LiteF;
 use Test::More;
+use Test::PDL;
 BEGIN {
 diag "ENV $_ = ", explain $ENV{$_}
   for qw(LD_LIBRARY_PATH DYLD_LIBRARY_PATH LDFLAGS CFLAGS CXXFLAGS LD_RUN_PATH
@@ -11,16 +12,6 @@ use PDL::Transform::Proj4;
 use Alien::proj;
 diag "Alien::proj version $Alien::proj::VERSION";
 
-sub tapprox
-{
-    my $x = shift;
-    my $y = shift;
-    my $d = abs($x - $y);
-    my $res = all($d < 1.0e-5);
-    diag "got:$x\nexpected:$y" if !$res;
-    $res;
-}
-
 my @version = eval { PDL::Transform::Proj4::proj_version() };
 is $@, '', 'proj_version no die';
 diag "PROJ version: @version";
@@ -38,10 +29,10 @@ my $xy_exp = double [
 ];
 
 my ($xy) = PDL::Transform::Proj4::fwd_transform($lonlat, $proj);
-ok( tapprox( $xy, $xy_exp ) );
+is_pdl $xy, $xy_exp;
 
 my ($lonlat2) = PDL::Transform::Proj4::inv_transform($xy, $proj);
-ok( tapprox( $lonlat2, $lonlat ) );
+is_pdl $lonlat2, $lonlat;
 
 # Do the corners of a cyl eq map, and see what we get...
 my $cyl_eq = "+proj=eqc +lon_0=0 +datum=WGS84";
@@ -54,7 +45,7 @@ my $xy3_exp = double [
 ];
 
 my ($xy3) = PDL::Transform::Proj4::fwd_transform($lonlat3, $cyl_eq);
-ok( tapprox( $xy3, $xy3_exp ) );
+is_pdl $xy3, $xy3_exp;
 
 $lonlat = ((xvals( double, 35, 17 ) - 17.0) * 10.0)->cat(
   (yvals( double, 35, 17 ) - 8.0) * 10.0)->mv(2,0);