#!/usr/bin/python

import sys
sys.path.append( "/home/server_b/modules/fahanalysis/" )
sys.path.append( "/home/server_b/modules/fahanalysis/FAHAnalysis3/Extras/" )
import FAHAnalysis3.AnalysisDriver as AnalysisDriver
import FAHAnalysis3.Helpers.Lock as Lock
from FAHAnalysis3.ToolKit import *
from FAHAnalysis3.ToolKit.BackupTool import PSCBackupTool
from FAHAnalysis3.Extras import TestCriteria
from FAHAnalysis3.Extras.TestCriteria import Definition 

sys.path.append( "/home/server_b/modules/fah_adaptive_sampling/adaptiveSampling/" )
from ConfigurationParser import ProjectConfigurationFile

myLock = Lock( ".analysis3060.lock" )

database = "project3060"
user = "server"
host = "vspmf92"
conffile = "/home/server_b/server2/densign/3060/project3060.conf" 
conf = ProjectConfigurationFile( conffile )

g_rms_bin = "/home/server_b/fahanalysis/binaries/g_rms"
g_gyrate_bin = "/home/server_b/fahanalysis/binaries/g_gyrate"
g_dist_bin = "/home/server_b/fahanalysis/binaries/g_dist"
g_mindist_bin = "/home/server_b/fahanalysis/binaries/g_mindist"

analyzer = AnalysisDriver.Analyzer( 
		datapath = conf.datapath,
		numruns = conf.nruns,
		numclones = conf.nclones,
		numsnapshots = conf.numdumps,
		sqlhost = host,
		sqluser = user,
		dbname = database, 
		filterxtc = True )

analyzer.dryrun = 0
analyzer.logsql = 0
analyzer.forcelastredo = 0
analyzer.testmode = 0
analyzer.trajext = "xtc"

# TOOLS
analyzer.addtool( AnalysisTimeTool() )

rg_tool = RGTool( indexfile="prot.ndx", binary = g_gyrate_bin, tprfile = "native.tpr" )
analyzer.addtool( rg_tool )

rmsdca_tool = RMSDTool33( "native.tpr", 0, indexfile = "calpha.ndx", binary = g_rms_bin )
rmsdca_tool.name = "RMSD-CA"
rmsdca_tool.fields = ( ("rmsdca", "float", "C alpha rmsd (Angstroms)"), )
analyzer.addtool( rmsdca_tool )

rmsdca_beta_turn = RMSDTool33( "native.tpr", 0, indexfile = "beta.ndx", binary = g_rms_bin )
rmsdca_beta_turn.name = "RMSD-CA of beta turn"
rmsdca_beta_turn.fields = ( ( "rmsdcaBT", "float", "C alpha rmsd of beta turn (Angstroms)" ), )
analyzer.addtool( rmsdca_beta_turn )

rmsdca_helix = RMSDTool33( "native.tpr", 0, indexfile = "alpha.ndx", binary = g_rms_bin )
rmsdca_helix.name = "RMSD-CA of helix"
rmsdca_helix.fields = ( ( "rmsdcaH", "float", "C alpha rmsd of helix (Angstroms)" ), )
analyzer.addtool( rmsdca_helix )

nat_contacts_tool = CANativeContactsTool( refStruct = "native.gro",
		tprFile = "native.tpr",
		cutoffDistance = 6.5,
		minGroupDist = 2,
		binary = "/home/server_b/fahanalysis/binaries/CAcontactMap" )
nat_contacts_tool.fields = (("natcon", "int", "number of native contacts") , )
analyzer.addtool( nat_contacts_tool )

tot_contacts_tool = CAContactsTool( tprFile = "native.tpr", 
		cutoffDistance = 6.5,
		minGroupDist = 2,
		binary = "/home/server_b/fahanalysis/binaries/CAcontactMap" )
tot_contacts_tool.fields = (("totcon", "int", "total number of contacts" ), )
analyzer.addtool( tot_contacts_tool )


# backup
backuptool = PSCBackupTool( "densign/project3060" )
backuptool.minchunksize = 0
backuptool.exclude_tpr = True
analyzer.addtool( backuptool )

dssp=DSSPTool( 23 ) 
analyzer.trjconv_indexgroup="0 0"
analyzer.trjconv_indexfile="index.ndx"
analyzer.addtool( dssp )

analyzer.addtool( PeriodicImageDistanceTool33(tpr = "native.tpr", indexgroup=0, indexfile="index.ndx", binary=g_mindist_bin ) )

#
"""
rclist=[]
for r in range( 11, conf.nruns ):
	for c in range( conf.nclones ):
		rclist.append( [r,c] )
"""
analyzer.go()

# folded criteria: from fitting histograms of these three things, and comparison with Young Min's data
crt1 = 3.15314  # rmsdca
crt2 = 0.994281 # rmsdcaBT
crt3 = 1.23099  # rmsdcaH
crt4 = 28.9 # native contacts (A)
crt5 = 33.7 # native contacts (B)

# show the criteria, and write to a file
FILE = open( "criteria", "a" )
summary = ( "rmsdca <= %3.3f angstroms\n" % crt1, "rmsdcaBT <= %3.3f angstroms\n" % crt2, "rmsdcaH <= %3.3f angstroms\n" % crt3)
FILE.writelines( summary )

# now add some definitional stuff
testCriteria = TestCriteria( dbname = database , sqlhost = host, sqluser = user )
# these correspond to the definition of 'folded' in the paper
definition1 = Definition( criteriaList = ("rmsdca <= %s" % crt1 , ) , updateCol = ("def1", 1 ) )
definition2 = Definition( criteriaList = ("rmsdcaBT <= %s" % crt2 , ) , updateCol = ("def2", 1 ) )
definition3 = Definition( criteriaList = ("rmsdcaH <= %s" % crt3 , ) , updateCol = ("def3", 1 ) )
definition4 = Definition( criteriaList = ("natcon > %s" % crt4, ) , updateCol = ("def4a", 1 ) )
definition5 = Definition( criteriaList = ("natcon > %s" % crt5, ) , updateCol = ("def4b", 1 ) )
 

testCriteria.addDefinition( definition1)# doRecip=True )
testCriteria.addDefinition( definition2)# doRecip=True )
testCriteria.addDefinition( definition3)# doRecip=True )
testCriteria.addDefinition( definition4)# doRecip=True )
testCriteria.addDefinition( definition5)# doRecip=True )

testCriteria.go()
myLock.unlock()