#!/usr/bin/awk -f
######################################################################################
# convert 3D 1H, 1HN, N NOE PCK file to con format
# usage: 3dhhn.awk null 3dnoe.PCK
######################################################################################

BEGIN { 
 skip = 1
 deleted = 0
 cutoff  = 0

 W_DIST = 4.0
 W_NEG  = 1.0 
 W_POS  = 0.9 

 M_DIST = 3.0
 M_NEG  = 1.0
 M_POS  = 0.9

 S_DIST = 2.5
 S_NEG  = 1.0
 S_POS  = 0.9

 WM     = 4.0
 MS     = 3.0

#MAX_INTENSITY = 1.03e07
MAX_INTENSITY = 9.57e05
MIN_INTENSITY = 1.38e04

CUTOFF_DIST = 4.8

#K = 0.34
K = (CUTOFF_DIST - 1.8)/(-1.0 + (MAX_INTENSITY/MIN_INTENSITY)^0.3333)
#print K
}

######################################################################################

######################################################################################
# Read in exclude defs definitions
######################################################################################
FILENAME == ARGV[1] { 
         if  ($1 == "#" || NF == 0) next
         exclude[$1] = 1
         next
}
######################################################################################
# skip over "header"
######################################################################################
FILENAME == ARGV[2] && skip { 
 	 while ($2 != "PkID")
		{
		 getline
		 }
	 getline
	 skip = 0
}

######################################################################################
function mk_ass(ass) {
 gsub(/[,;]/,".",ass)
 ns = split(ass, aa, ".")
 ASS = sprintf( "%s ", aa[3])
 return substr(aa[1],2,length(aa[1]))
}

function mk_ass2(ass) { 
 gsub(/[,;]/,".",ass)
 ns = split(ass, aa, ".")
 ASS = sprintf( "%s ", aa[2])
 return substr(aa[1],2,length(aa[1]))
}
######################################################################################
function pr_output( swap, res1, ass1, res2, ass2, dist, neg, pos) { 
 if (swap)
 	printf "assign (resid %-3s and name %-4s) (resid %-3s and name %-4s) ",res2, ass2, res1, ass1

 else
	printf "assign (resid %-3s and name %-4s) (resid %-3s and name %-4s) ",res1, ass1, res2, ass2

# output constraints
 printf "%4.1f %4.1f %4.1f\n", dist, neg, pos

}
######################################################################################
{
# screen for unwanted peak ids

 if ($1 in exclude)
        {if (deleted == 0)
	 	print $0 > "noe.deleted"
	 else
	 	print $0 >> "noe.deleted"
	 deleted++
         next
        }

#skip unassigned
 if (match($7,  /[ACDEFGHIKLMNPQRSTVWY][0-9]*\./) == 0 || \
     match($8,  /[ACDEFGHIKLMNPQRSTVWY][0-9]*\./) == 0 )
        {next}

 res1 = mk_ass($7)
 ass1 = ASS
 res2 = mk_ass2($8)
 ass2 = ASS

# printf ">>> res1=%s, ass1=%s, res2=%s, ass2=%s\n",res1,ass1,res2,ass2

# always put smallest residue id first, or "lower" proton identifier first for
# intra-residue
 swap = (int(res1) > int(res2)) || (int(res1) == int(res2) && ass1 > ass2)

# calc constraint limits and print output
 
 absi = ($6 < 0.0) ? -1.0*$6 : $6
 absi = MAX_INTENSITY / absi
 absi = (absi)^0.3333
 absi = 1.8 + (-1.0 + absi) * K

if(res1==res2 && ass1==ass2) {next}

 if (absi > CUTOFF_DIST)
	{if (cutoff == 0)
		{print $0 > "noe.cutoff"}
	 else
		{print $0 >> "noe.cutoff"}
	 cutoff++
	}
 else if (absi <= CUTOFF_DIST && absi >= WM)
        {pr_output(swap, res1, ass1, res2, ass2, absi, 0.40*absi, 0.40*absi)}
 else if (absi < WM && absi >= MS)
        {pr_output(swap, res1, ass1, res2, ass2, absi, 0.40*absi, 0.40*absi)}
 else if (absi < MS && absi >= 2.2)
        {pr_output(swap, res1, ass1, res2, ass2, absi, 0.40*absi, 0.40*absi)}
 else if (absi < 2.2)
         {pr_output(swap, res1, ass1, res2, ass2, absi, 0.40*absi, 0.40*absi)}
}