我写的用来计算一个二聚体蛋白之间相靠近的tk/tcl的script
有什么批评和建议请尽管提出来!:)
proc Inter_Dist_Cal {} {
set di "1Q2W.pdb"
puts -nonewline "Enter PDB file (press \"d\" as \"$di\"): "
flush stdout
set infile [gets stdin]
if {[string equal -nocase $infile "d"]} {set infile $di}
set PDBname [lindex [split $infile .] 0]
set do1 "($PDBname)All_Inter_Con.txt"
puts "All Inter Contact will be \"$do1\" "
flush stdout
set outfile1 $do1
set do2 "($PDBname)Res_Inter_Con.txt"
puts "Residual Inter Contact will be \"$do2\" "
flush stdout
set outfile2 $do2
set do3 "($PDBname)Num_Res_Inter_Con.txt"
puts "The number of Inter Contact for each Residue will be \"$do3\" "
flush stdout
set outfile3 $do3
set dd "10"
puts -nonewline "Enter max distance (press \"d\" as \"$dd\"): "
flush stdout
set indist [gets stdin]
if {[string equal -nocase $indist "d"]} {set indist $dd}
set in [open $infile r]
set n1 0
set n2 0
while {[gets $in line]>=0} {
#set line [lindex [split $line "\""] 0]
if {([lindex $line 0]==" ")||($line=="")||([lindex [split [lindex $line 0] ""] 0]=="#")}
if {[lindex $line 0]!="ATOM"}
if {[lindex $line 4]=="A"} {
set 1id($n1) [lindex $line 1]
set 1atom($n1) [lindex $line 2]
set 1resty($n1) [lindex $line 3]
set 1resno($n1) [lindex $line 5]
set 1x($n1) [lindex $line 6]
set 1y($n1) [lindex $line 7]
set 1z($n1) [lindex $line 8]
incr n1
}
if {[lindex $line 4]=="B"} {
set 2id($n2) [lindex $line 1]
set 2atom($n2) [lindex $line 2]
set 2resty($n2) [lindex $line 3]
set 2resno($n2) [lindex $line 5]
set 2x($n2) [lindex $line 6]
set 2y($n2) [lindex $line 7]
set 2z($n2) [lindex $line 8]
incr n2
}
}
close $in
set out1 [open $outfile1 w]
puts $out1 "# [clock format [clock seconds]] BY [lindex [exec who] end-4]"
puts $out1 ""
puts $out1 [format "# %4s %5s %5s - %6s %5s %5s %9s" Mol1 Res Atom Mol2 Res Atom Distance]
puts $out1 ""
for {set i 0} {incr i} {
for {set j 0} {incr j} {
if {[set dis [expr pow(pow(x($i)-x($j),2)+pow(y($i)-y($j),
2)+pow(z($i)-z($j),2),0.5)]]<=$indist} {
set fm [format "%6d %5s %5s - %6d %5s %5s %9.3f" resno($i) resty($i) atom($i)
resno($j) resty($j) atom($j) $dis]
puts $out1 $fm
set temp1 resno($i)
if {[info exists ACon($temp1)]} {
incr ACon($temp1)
} else {set ACon($temp1) 1}
set temp2 [expr resno($i)*1000+resno($j)]
if {[info exists ABCon($temp2)]} {
incr ABCon($temp2)
} else {set ABCon($temp2) 1}
}
}
}
close $out1
set out2 [open $outfile2 w]
puts $out2 "# [clock format [clock seconds]] BY [lindex [exec who] end-4]"
puts $out2 ""
puts $out2 [format "# %4s %5s %6s %5s %9s" Mol1 Res Mol2 Res Contact]
puts $out2 ""
for {set i 0} {incr i} {
for {set j 0} {incr j} {
if {[info exists ABCon([expr resno($i)*1000+resno($j)])]} {
set fm [format "%6d %5s %6d %5s %9d" resno($i) resty($i) resno($j) resty($j)
$ABCon([expr resno($i)*1000+resno($j)])]
puts $out2 $fm
unset ABCon([expr resno($i)*1000+resno($j)])
}
}
}
close $out2
set out3 [open $outfile3 w]
puts $out3 "# [clock format [clock seconds]] BY [lindex [exec who] end-4]"
puts $out3 ""
puts $out3 [format "# %4s %5s %9s" Mol1 Res Contact]
puts $out3 ""
for {set i 0} {incr i} {
if {[info exists ACon(resno($i))]} {
set fm [format "%6d %5s %9d" resno($i) resty($i) $ACon(resno($i))]
puts $out3 $fm
unset ACon(resno($i))
}
}
close $out3
}