#
# Read a SIMWALK2 format haplotype file, return a digraph in the dot language 
# call dot to make pedigree drawing
# 
AWK=gawk
DOTEXE=dot
if [ "$1" = "" ] || [ "$1" = "-h" ] || [ "$1" = "--help" ]
then
        echo "Usage: drawhap.sh [--output <format>] <SIMWALK2 haplotype file> ..."
        exit
fi

OUTFORM=ps
switch=1
while [ $switch = 1 ]
do case "$1" in
  -o) OUTFORM=$2; shift 2 ;;
  --output) OUTFORM=$2; shift 2 ;;
  -*) echo "Flag not recognised: $1" ; shift ;;
  *) switch=0 ;;
  esac
done


for fil in $*
do
  ped=`$AWK '/The results for the pedigree named:/ {print $NF}' $fil`
  echo "Pedigree $ped"
  $AWK -v ped=$ped '
  BEGIN { shape["M"]="box,regular=1"
          shape["F"]="circle"
          shade["AFFECTED"]="grey"
          shade["NORMAL"]="white"
          shade["--"]="white"
          palette[1]="yellow1"
          palette[2]="turquoise"
          palette[3]="tan"
          palette[4]="lightblue1"
          palette[5]="wheat"
          palette[6]="thistle2"
          palette[7]="yellow2"
          palette[8]="lightblue3"
          palette[9]="wheat3"
          palette[10]="violet"
          palette[11]="steelblue"
          palette[12]="grey"
          for(i=1;i<=100;i++) {
            colour[i,"-"]="transparent"
          }
          for(i=1;i<=14;i++) getline
  }
  substr($0,1,1)=="M"  { 
          id=$2
          sex[id]=$5
          aff[id]=$6
          if($3!="--") {
            kids[$3]++
            kids[$4]++
            marriage[$3,$4]++
            child[$3,$4,marriage[$3,$4]]=id
          }
          gsub("@","-")
          nloci=NF-6
          blanks=0
          for(i=7;i<=NF;i++) if ($i=="-") blanks++
          if (nloci>blanks) {
            tbl="<table border=\"0\" cellborder=\"0\">\n"
            for(i=7;i<=NF;i++) {
              loc=i-6
              if (colour[loc, $i]=="") {
                allele[loc]++
                colour[loc, $i]=palette[allele[loc]]
              }
              tbl=tbl "<tr><td bgcolor=\"" colour[loc,$i] "\">" $i "</td></tr>\n"
            }
            tbl=tbl "</table>"
          }else{
            tbl=""
          }
          hap1[id]=tbl
          getline
          split($2,rec)
          getline
          gsub("@","-")
          blanks=0
          for(i=7;i<=NF;i++) if ($i=="-") blanks++
          if (nloci>blanks) {
            tbl="<table border=\"0\" cellborder=\"0\">\n"
            for(i=7;i<=NF;i++) {
              loc=i-6
              if (colour[loc, $i]=="") {
                allele[loc]++
                colour[loc, $i]=palette[allele[loc]]
              }
              tbl=tbl "<tr><td bgcolor=\"" colour[loc,$i] "\">" $i "</td></tr>\n"
            }
            tbl=tbl "</table>"
          }else{
            tbl=""
          }
          hap2[id]=tbl
  }
  END   { print "digraph Ped_" ped " {"
          print "# page =\"8.2677165,11.692913\" ;"
          print "size =\"8.2677165,11.692913\" ;"
          print "# ratio =\"auto\" ;"
          print "mincross = 2.0 ;"
          print "# label=\"Pedigree " ped "\" ;"
          print "# rotate=90 ;"
          print "ranksep=1.5 ;"
          for(s in sex) if (!(s in hap1)) {
            print "\"" s "\" [shape=" shape[sex[s]] ","  \
                  " style=filled,fillcolor=" shade[aff[s]] "] ;"
          }
          for(h in hap1) {
            if (kids[h]>0) {
              print "subgraph " h " { "
              print "  \"" h "\" [shape=" shape[sex[h]] ","  \
                    "   style=filled,fillcolor=" shade[aff[h]] "] ;"
              print "  \"" h "hap1\" [rank=same, shape=plaintext \\"   
              print "           label=<" hap1[h] ">, margin=0.01] ;"
              print "  \"" h "hap2\" [rank=same, shape=plaintext \\"   
              print "           label=<" hap1[h] ">] ;"
              print "  \"" h "mid\" [rank=same, shape=circle, label=\"\", height=0, width=0]"   
              print "  \"" h "\" -> \"" h "hap1\" [color=transparent, weight=500] ;"
              print "  \"" h "\" -> \"" h "mid\" [dir=none, weight=1000] ;"
              print "  \"" h "\" -> \"" h "hap2\" [color=transparent, weight=500] ;"
              print "  \"" h "hap1\" -> \"" h "mid\" [color=transparent, weight=100] ;"
              print "  \"" h "hap2\" -> \"" h "mid\" [color=transparent, weight=100] ;"
              print "  \"" h "base\" [shape=circle, label=\"\", height=0, width=0]"   
              print "  \"" h "mid\" -> \"" h "base\" [dir=none, weight=1000] ;"
              print "}"
            }else{
              print "subgraph {"
              print "  \"" h "\" [shape=" shape[sex[h]] ","  \
                    "   style=filled,fillcolor=" shade[aff[h]] "] ;"
              print "  \"" h "\" -> \"" h "hap1\" [color=transparent, weight=100] ;"
              print "  \"" h "\" -> \"" h "hap2\" [color=transparent, weight=100] ;"
              print "\"" h "hap1\" [rank=same, shape=plaintext \\"   
              print "         label=<" hap1[h] ">] ;"
              print "\"" h "hap2\" [rank=same, shape=plaintext \\"   
              print "         label=<" hap1[h] ">] ;"
              print "}"
            }
          }
          for(m in marriage) {
            n=split(m,par,"\034")
            mating="\"" par[1] "x" par[2] "\""
            print mating "[shape=diamond,style=filled," \
                  "label=\"\",height=.1,width=.1] ;"
            print "\"" par[1] "base\" -> " mating " [dir=none, weight=1] ;"
            print "\"" par[2] "base\" -> " mating " [dir=none, weight=1] ;"
            for(k=1;k<=marriage[par[1],par[2]];k++) {
              print  mating " -> \"" child[par[1],par[2],k] "\"" \
                     " [dir=none, weight=2] ;"
            }
          }
          print "}"
  }' $fil > $ped.dot 
  $DOTEXE -T$OUTFORM $ped.dot -o $ped.$OUTFORM
done

