#!/usr/bin/perl -T


# jigsaw puzzle assembler, for e.g. Google Map tiles
# A.Daviel, TRIUMF, Aug 2007

# assembles square image tiles into a Fig drawing by matching up the edges.
# Edge matching is done by computing sum of distances in RGB colorspace between
# corresponding pixels in 1-pixel-wide edge strips.

# assumes current directory contains tile JPEGs. Tiles must currently be square 
# and must be all the same way up and the same size
# A regex pattern ($filepat) is used to identify the correct tiles

# a number of temporary files (*.rgb) are created and left around. If present, subsequent
# script execution is faster

# Note - Google maps and satellite images are copyright. Depending on jurisdiction,
# you may be able to make single copies for private use or research purposes


$T = 6000 ; # threshold for match - distance should be less than this for matches and more for mismatch
$size = 256 ; # size of square tiles in pixels
$filepat = '^(tqt.................jpg)$' ; # pattern to match tile name
#             tqtrtrtttstrqqqsrsr.jpg
$filepat = '^(tq[\w]{17}.jpg)$' ; # pattern to match tile name

$dbg = shift ; # debug flag

$size3 = $size * 3 ;
$sizem1 = $size - 1 ;
$sizex1 = $size . 'x1' ; 

$ENV{PATH} = '/usr/bin' ; # need "convert" from ImageMagick package

opendir (DIR1,".") or die "Opendir failed"  ;
foreach $_ (readdir (DIR1)) {
  if ($dbg>1) { print "$_\n" ; }
  unless (/$filepat/) { next ; }     # pattern to match tile name
  $file1 = $1 ;
  push (@file,$file1) ;
  unless (-f "$file1.L.rgb") {
    system("convert -crop 1x$size+0+0 $file1 $file1.L.rgb") ;
  }
  unless (-f "$file1.R.rgb") {
    system("convert -crop 1x$size+$sizem1+0 $file1 $file1.R.rgb") ; 
  }
  unless (-f "$file1.T.rgb") {
    system("convert -crop $sizex1+0+0 $file1 $file1.T.rgb") ; 
  }
  unless (-f "$file1.B.rgb") {
    system("convert -crop $sizex1+0+$sizem1 $file1 $file1.B.rgb") ; 
  }
  open (IN,"$file1.R.rgb") or die "Need $file1.R.rgb" ;
  $n = read (IN,$AR,$size3) ;
  unless ($n == $size3) { die "AR not $size3 bytes\n" ; }
  open (IN,"$file1.L.rgb") or die "Need $file1.L.rgb" ;
  $n = read (IN,$AL,$size3) ;
  unless ($n == $size3) { die "AL not $size3 bytes\n" ; }
  open (IN,"$file1.T.rgb") or die "Need $file1.T.rgb" ;
  $n = read (IN,$AT,$size3) ;
  unless ($n == $size3) { die "AT not $size3 bytes\n" ; }
  open (IN,"$file1.B.rgb") or die "Need $file1.B.rgb" ;
  $n = read (IN,$AB,$size3) ;
  unless ($n == $size3) { die "AB not $size3 bytes\n" ; }

  $D13m = $D24m = $D57m = $D68m = 999999;
  opendir (DIR,".") or die "Opendir failed"  ;
  foreach $_ (readdir (DIR)) {
    unless (/$filepat/) { next ; }     # pattern to match tile name
    $file = $1 ;
    if ($file eq $file1) { next ; }

    unless (-f "$file.L.rgb") {
      if ($dbg>1) { print "convert -crop 1x$size+0+0 $file $file.L.rgb\n" ; }
      system("convert -crop 1x$size+0+0 $file $file.L.rgb") ;
    }
    unless (-f "$file.R.rgb") {
      if ($dbg>1) { print "convert -crop 1x$size+$sizem1+0 $file $file.R.rgb\n" ; }
      system("convert -crop 1x$size+$sizem1+0 $file $file.R.rgb") ;
    }
    unless (-f "$file.T.rgb") {
      if ($dbg>1) { print "convert -crop $sizex1+0+0 $file $file.T.rgb\n" ; }
      system("convert -crop $sizex1+0+0 $file $file.T.rgb") ;
    }
    unless (-f "$file.B.rgb") {
      if ($dbg>1) { print "convert -crop $sizex1+0+$sizem1 $file $file.B.rgb\n" ; }
      system("convert -crop $sizex1+0+$sizem1 $file $file.B.rgb") ;
    }

    open (IN,"$file.L.rgb") or die "Need $file.L.rgb" ;
    $n = read (IN,$BL,$size3) ;
    unless ($n == $size3) { print "$file BL not $size3 bytes\n" ; next ; }
    open (IN,"$file.R.rgb") or die "Need $file.R.rgb" ;
    $n = read (IN,$BR,$size3) ;
    unless ($n == $size3) { print "$file BR not $size3 bytes\n" ; next ; }
    open (IN,"$file.T.rgb") or die "Need $file.T.rgb"  ;
    $n = read (IN,$BT,$size3) ;
    unless ($n == $size3) { print "$file BT not $size3 bytes\n" ; next ; }
    open (IN,"$file.B.rgb") or die "Need $file.B.rgb" ;
    $n = read (IN,$BB,$size3) ;
    unless ($n == $size3) { print "$file BB not $size3 bytes\n" ; next ; }
    @pix1 = split(//,$AR) ; @pix2 = split(//,$AL) ; @pix5 = split(//,$AT) ; @pix6 = split(//,$AB) ;
    @pix3 = split(//,$BL) ; @pix4 = split(//,$BR) ; @pix7 = split(//,$BB) ; @pix8 = split(//,$BT) ;

    $D13 = $D24 = $D57 = $D68 = 0 ;
    for ($i=0;$i<$size;$i++) {
      $r1 = ord(shift(@pix1)) ; $g1 = ord(shift(@pix1)) ; $b1 = ord(shift(@pix1)) ;
      $r2 = ord(shift(@pix2)) ; $g2 = ord(shift(@pix2)) ; $b2 = ord(shift(@pix2)) ;
      $r3 = ord(shift(@pix3)) ; $g3 = ord(shift(@pix3)) ; $b3 = ord(shift(@pix3)) ;
      $r4 = ord(shift(@pix4)) ; $g4 = ord(shift(@pix4)) ; $b4 = ord(shift(@pix4)) ;
      $r5 = ord(shift(@pix5)) ; $g5 = ord(shift(@pix5)) ; $b5 = ord(shift(@pix5)) ;
      $r6 = ord(shift(@pix6)) ; $g6 = ord(shift(@pix6)) ; $b6 = ord(shift(@pix6)) ;
      $r7 = ord(shift(@pix7)) ; $g7 = ord(shift(@pix7)) ; $b7 = ord(shift(@pix7)) ;
      $r8 = ord(shift(@pix8)) ; $g8 = ord(shift(@pix8)) ; $b8 = ord(shift(@pix8)) ;

      $d13 = sqrt( ($r1-$r3)**2 + ($b1-$b3)**2 + ($g1-$g3)**2 ) ;
      $D13+= $d13 ;
      $d24 = sqrt( ($r2-$r4)**2 + ($b2-$b4)**2 + ($g2-$g4)**2 ) ;
      $D24+= $d24 ;

      $d57 = sqrt( ($r5-$r7)**2 + ($b5-$b7)**2 + ($g5-$g7)**2 ) ;
      $D57+= $d57 ;
      $d68 = sqrt( ($r6-$r8)**2 + ($b6-$b8)**2 + ($g6-$g8)**2 ) ;
      $D68+= $d68 ;

    }
    $D13 = int($D13) ; $D24 = int($D24) ; $D57 = int($D57) ; $D68 = int($D68) ;
    if ($dbg) {
      print "$file AR-BL $D13 AL-BR $D24 AT-BB $D57 AB-BT $D68\n" ;
    }
#    unless(defined($D13m)) { $D13m = $D13 ; }
#    unless(defined($D24m)) { $D24m = $D24 ; }
#    unless(defined($D57m)) { $D57m = $D57 ; }
#    unless(defined($D68m)) { $D68m = $D68 ; }
    if ($D13<$D13m) { $D13m = $D13 ; $left = $file ; }
    if ($D24<$D24m) { $D24m = $D24 ; $right= $file ; }
    if ($D57<$D57m) { $D57m = $D57 ; $botm = $file ; }
    if ($D68<$D68m) { $D68m = $D68 ; $top  = $file ; }
  }
  print "$file1 is\n" ; 
  if ($D13m < $T) { $leftof{$file1} = $left ;   print "  Left of  $left ($D13m)\n" ; }
  if ($D24m < $T) { $rightof{$file1} = $right ; print "  Right of $right ($D24m)\n" ; }
  if ($D68m < $T) { $above{$file1} = $top ;     print "  Above    $top ($D68m)\n" ; }
  if ($D57m < $T) { $below{$file1} = $botm ;    print "  Below    $botm ($D57m)\n" ;}
}
$x{$file[0]} = 10 ; $y{$file[0]} = 10 ; 

for ($k=0;$k<4;$k++) {
  foreach $file (@file) {
    if ($x{$file}) {
      if ($leftof{$file}) {
	$x{$leftof{$file}} = $x{$file} + 1 ;
	$y{$leftof{$file}} = $y{$file} ;
      }
      if ($rightof{$file}) {
	$x{$rightof{$file}} = $x{$file} - 1 ;
	$y{$rightof{$file}} = $y{$file} ;
      }
      if ($above{$file}) {
	$x{$above{$file}}  = $x{$file} ;
	$y{$above{$file}}  = $y{$file} +1 ;
      }
      if ($below{$file}) {
	$x{$below{$file}}  = $x{$file} ;
	$y{$below{$file}} = $y{$file}  - 1 ;
      }
    }
  }
}

open (OUT,">tiles.fig") ;
print OUT<<EOT;
#FIG 3.2
Landscape
Center
Inches
Letter
100.00
Single
-2
# created by $0
1200 2
EOT
$scale = 3840 ;
foreach $file (@file) {
  $x1 = $x{$file} * $scale ; $y1 = $y{$file} * $scale ;
  $x2 = $x1 + $scale ; $y2 = $y1 + $scale ;
  print OUT<<EOT;
# $file
2 5 0 1 0 -1 60 -1 -1 0.000 0 0 -1 0 0 5
	0 $file
	$x1 $y1 $x2 $y1 $x2 $y2 $x1 $y2 $x1 $y1
EOT
}