-
Notifications
You must be signed in to change notification settings - Fork 2
/
rmsdtt.tcl
1803 lines (1560 loc) · 55.7 KB
/
rmsdtt.tcl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
#
# RMSD Trajectory Tool v4.0
#
# A GUI interface for RMSD alignment and analysis of trajectories
#
# http://physiology.med.cornell.edu/faculty/hweinstein/vmdplugins/rmsdtt
# Author
# ------
# Luis Gracia, PhD
# Weill Medical College, Cornel University, NY
# lug2002@med.cornell.edu
# This plugin is based on the RMSD Tool plugin written by
# Alexander Balaeff, PhD, TB Group, UIUC
# Pascal Mercier, PhD, Biochemistry Dept, Univ. of Alberta, Edmonton, Canada
# John Stone, TB Group, UIUC
# Documentation
# -------------
# See the index.html file distributed with this file. For more update documentation see
# http://physiology.med.cornell.edu/faculty/hweinstein/vmdplugins/rmsdtt
package provide rmsdtt 4.0
namespace eval ::rmsdtt:: {
namespace export rmsdtt
}
proc rmsdtt_tk_cb {} {
# This gets called by VMD the first time the menu is opened.
#variable foobar
# Don't destroy the main window, because we want to register the window
# with VMD and keep reusing it. The window gets iconified instead of
# destroyed when closed for any reason.
#set foobar [catch {destroy $::rmsdtt::w }] ;# destroy any old windows
# start RMSDTT
::rmsdtt::rmsdtt
return $rmsdtt::w
}
proc rmsdtt::rmsdtt {} {
variable w
variable bb_only 0
variable trace_only 0
variable noh 1
variable swap_use 1
variable swap_sw 0
variable swap_type all
variable swap_print 1
variable rmsd_base top
variable traj_sw 1
variable traj_ref 0
variable traj_all 0
variable skip_sw 0
variable skip_start 0
variable skip_end end
variable skip_steps 1
variable time_sw 0
variable time_ini 0.0
variable time_step 1.0
variable save_sw 0
variable save_file trajrmsd.dat
variable plot_sw 0
variable plot_program multiplot
variable colorize 1
variable bb_def "C CA N"
variable stats 1
variable datalist
variable datatot
variable equiv_sw 0
variable equiv_byres 0
variable equiv_atom CA
variable equiv_lines {}
variable equiv_mol top
variable weighted_sw 0
variable weighted_field user
variable weighted_mol ref
# If already initialized, just turn on
if { [winfo exists .rmsdtt] } {
wm deiconify $w
return
}
# GUI look
#option add *rmsdtt.*font {Helvetica 9}
#option add *rmsdtt.top.left.sel.background white
#option add *rmsdtt.*Text.background white
option add *rmsdtt.*borderWidth 1
option add *rmsdtt.*Button.padY 0
option add *rmsdtt.*Menubutton.padY 0
# Main window
set w [toplevel ".rmsdtt"]
wm title $w "WMC PhysBio - RMSD Trajectory Tool"
# Menu
frame $w.menubar -relief raised -bd 2
pack $w.menubar -fill x
menubutton $w.menubar.file -text "File" -menu $w.menubar.file.menu -underline 0 -pady 2
menu $w.menubar.file.menu -tearoff no
$w.menubar.file.menu add command -label "Save data..." -command "[namespace current]::SaveDataBrowse data" -underline 0
$w.menubar.file.menu add command -label "Save summary..." -command "[namespace current]::SaveDataBrowse summary" -underline 0
$w.menubar.file.menu add command -label "Plot data" -command [namespace current]::doPlot -underline 0
pack $w.menubar.file -side left
menubutton $w.menubar.options -text "Options" -menu $w.menubar.options.menu -underline 0 -pady 2
menu $w.menubar.options.menu -tearoff no
$w.menubar.options.menu add cascade -label "Backbone def..." -menu $w.menubar.options.menu.bbdef -underline 0
menu $w.menubar.options.menu.bbdef
$w.menubar.options.menu.bbdef add radiobutton -label "C CA N" -variable [namespace current]::bb_def -value "C CA N"
$w.menubar.options.menu.bbdef add radiobutton -label "C CA N O" -variable [namespace current]::bb_def -value "C CA N O"
$w.menubar.options.menu add cascade -label "Plotting program..." -menu $w.menubar.options.menu.plot -underline 0
menu $w.menubar.options.menu.plot -tearoff no
$w.menubar.options.menu.plot add radiobutton -label "Multiplot (all)" -variable [namespace current]::plot_program -value "multiplot" -underline 0
$w.menubar.options.menu.plot add radiobutton -label "Xmgrace (Unix)" -variable [namespace current]::plot_program -value "xmgrace" -underline 0
$w.menubar.options.menu.plot add radiobutton -label "MS Excel (Windows)" -variable [namespace current]::plot_program -value "excel" -underline 4
$w.menubar.options.menu add checkbutton -label "Statistics" -variable [namespace current]::stats -underline 0
$w.menubar.options.menu add checkbutton -label "Colorize table" -variable [namespace current]::colorize -underline 0
pack $w.menubar.options -side left
menubutton $w.menubar.help -text "Help" -menu $w.menubar.help.menu -underline 0 -pady 2
menu $w.menubar.help.menu -tearoff no
$w.menubar.help.menu add command -label "About" -command [namespace current]::help_about
$w.menubar.help.menu add command -label "Help..." -command "vmd_open_url http://physiology.med.cornell.edu/faculty/hweinstein/vmdplugins/rmsdtt/index.html"
pack $w.menubar.help -side right
# top frame
frame $w.top
pack $w.top -side top -fill x
# Selection
frame $w.top.left -relief ridge
pack $w.top.left -side left -fill both -expand yes
text $w.top.left.sel -height 3 -width 25 -highlightthickness 0 -selectborderwidth 0 -exportselection yes -wrap word -relief sunken -bd 1
pack $w.top.left.sel -side top -fill both -expand yes
$w.top.left.sel insert end "protein"
bind $w.top.left.sel <Return> "[namespace current]::draw_equiv"
labelframe $w.top.left.mods -text "Selection Modifiers" -relief ridge -bd 2
pack $w.top.left.mods -side top -fill x
checkbutton $w.top.left.mods.bb -text "Backbone" -variable [namespace current]::bb_only -command "[namespace current]::ctrlbb bb"
checkbutton $w.top.left.mods.tr -text "Trace" -variable [namespace current]::trace_only -command "[namespace current]::ctrlbb trace"
checkbutton $w.top.left.mods.noh -text "noh" -variable [namespace current]::noh -command "[namespace current]::ctrlbb noh"
menubutton $w.top.left.mods.selectionhistory -menu $w.top.left.mods.selectionhistory.m -text "History" -relief raised -direction flush
menu $w.top.left.mods.selectionhistory.m -tearoff no
pack $w.top.left.mods.bb $w.top.left.mods.tr $w.top.left.mods.noh -side left -anchor w
pack $w.top.left.mods.selectionhistory -side right
# Draw lines for equivalent atoms
labelframe $w.top.left.equiv -text "Equivalent Atoms" -relief ridge -bd 2
pack $w.top.left.equiv -side top -fill x
checkbutton $w.top.left.equiv.0 -text "On/Off" -variable [namespace current]::equiv_sw -command [namespace current]::draw_equiv
checkbutton $w.top.left.equiv.byres -text "byres" -variable [namespace current]::equiv_byres -command [namespace current]::draw_equiv
label $w.top.left.equiv.atomlabel -text "atom:"
entry $w.top.left.equiv.atom -width 4 -textvariable [namespace current]::equiv_atom
pack $w.top.left.equiv.0 $w.top.left.equiv.byres $w.top.left.equiv.atomlabel $w.top.left.equiv.atom -side left -anchor w
bind $w.top.left.equiv.atom <Return> "[namespace current]::draw_equiv"
# Swap
if {[catch {package require swap} msg]} {
set swap_use 0
}
if {$swap_use} {
labelframe $w.top.left.swap -text "Swap Atoms" -relief ridge -bd 2
pack $w.top.left.swap -side top -fill x
checkbutton $w.top.left.swap.0 -text "On/Off" -variable [namespace current]::swap_sw -command [namespace current]::ctrlgui
menubutton $w.top.left.swap.type -relief raised -direction flush -textvariable [namespace current]::swap_type -menu $w.top.left.swap.type.menu
menu $w.top.left.swap.type.menu -tearoff no
checkbutton $w.top.left.swap.print -text "verbose" -variable [namespace current]::swap_print
button $w.top.left.swap.list -relief raised -text "List" -command [namespace code {::swap::list $swap_type}]
pack $w.top.left.swap.0 $w.top.left.swap.type $w.top.left.swap.print -side left -anchor w
pack $w.top.left.swap.list -side right
}
# Weighted rmsd
labelframe $w.top.left.weighted -text "Weights" -relief ridge -bd 2
pack $w.top.left.weighted -side top -fill x
checkbutton $w.top.left.weighted.0 -text "On/Off" -variable [namespace current]::weighted_sw -command [namespace current]::ctrlgui
label $w.top.left.weighted.mollabel -text "Mol:"
menubutton $w.top.left.weighted.mol -relief raised -direction flush -textvariable [namespace current]::weighted_mol -menu $w.top.left.weighted.mol.menu
menu $w.top.left.weighted.mol.menu -tearoff no
label $w.top.left.weighted.fieldlabel -text "Field:"
menubutton $w.top.left.weighted.field -relief raised -direction flush -textvariable [namespace current]::weighted_field -menu $w.top.left.weighted.field.menu
menu $w.top.left.weighted.field.menu -tearoff no
pack $w.top.left.weighted.0 $w.top.left.weighted.mollabel $w.top.left.weighted.mol $w.top.left.weighted.fieldlabel $w.top.left.weighted.field -side left -anchor w
foreach field [list user mass charge beta occupancy] {
$w.top.left.weighted.field.menu add radiobutton -value $field -label $field -variable [namespace current]::weighted_field
}
# Buttons
frame $w.top.right
pack $w.top.right -side right
frame $w.top.right.pushfr -relief ridge -bd 2
pack $w.top.right.pushfr -side top -fill x
button $w.top.right.pushfr.rmsd -text "RMSD" -relief raised -command [namespace current]::doRmsd -pady 2 -bd 2
button $w.top.right.pushfr.align -text "ALIGN" -relief raised -command [namespace current]::doAlign -pady 2 -bd 2
pack $w.top.right.pushfr.rmsd $w.top.right.pushfr.align -side left -fill x -expand yes
# Ref mol
labelframe $w.top.right.ref -text "Reference mol" -relief ridge -bd 2
pack $w.top.right.ref -side top -fill x
radiobutton $w.top.right.ref.0 -text "Top" -variable [namespace current]::rmsd_base -value "top" -command [namespace current]::ctrlgui
radiobutton $w.top.right.ref.1 -text "Average" -variable [namespace current]::rmsd_base -value "ave" -command [namespace current]::ctrlgui
radiobutton $w.top.right.ref.2 -text "Selected" -variable [namespace current]::rmsd_base -value "selected" -command [namespace current]::ctrlgui
pack $w.top.right.ref.0 $w.top.right.ref.1 $w.top.right.ref.2 -side left -expand yes
# Trajectory
labelframe $w.top.right.traj -text "Trajectory" -relief ridge -bd 2
pack $w.top.right.traj -side top
frame $w.top.right.traj.frames -relief ridge -bd 0
pack $w.top.right.traj.frames -side top -fill x
checkbutton $w.top.right.traj.frames.0 -text "On/Off" -variable [namespace current]::traj_sw -command [namespace current]::ctrlgui
label $w.top.right.traj.frames.reflabel -text "Frame ref:"
entry $w.top.right.traj.frames.ref -width 5 -textvariable [namespace current]::traj_ref
checkbutton $w.top.right.traj.frames.all -text "All" -variable [namespace current]::traj_all -command [namespace current]::ctrlgui
pack $w.top.right.traj.frames.0 $w.top.right.traj.frames.reflabel $w.top.right.traj.frames.ref -side left -anchor w -fill x
pack $w.top.right.traj.frames.all -side right
frame $w.top.right.traj.skip -relief ridge -bd 0
pack $w.top.right.traj.skip -side top -fill x
checkbutton $w.top.right.traj.skip.0 -text "Skip" -variable [namespace current]::skip_sw -command [namespace current]::ctrlgui
label $w.top.right.traj.skip.inilabel -text "Start:"
entry $w.top.right.traj.skip.ini -width 5 -textvariable [namespace current]::skip_start
label $w.top.right.traj.skip.endlabel -text "End:"
entry $w.top.right.traj.skip.end -width 5 -textvariable [namespace current]::skip_end
label $w.top.right.traj.skip.stepslabel -text "Steps:"
entry $w.top.right.traj.skip.steps -width 5 -textvariable [namespace current]::skip_steps
pack $w.top.right.traj.skip.0 $w.top.right.traj.skip.inilabel $w.top.right.traj.skip.ini $w.top.right.traj.skip.stepslabel $w.top.right.traj.skip.steps -side left -anchor w
# $w.top.right.traj.skip.endlabel $w.top.right.traj.skip.end
frame $w.top.right.traj.time -relief ridge -bd 0
pack $w.top.right.traj.time -side top -fill x
checkbutton $w.top.right.traj.time.0 -text "Time" -variable [namespace current]::time_sw -command [namespace current]::ctrlgui
label $w.top.right.traj.time.inilabel -text "Start:"
entry $w.top.right.traj.time.inival -width 6 -textvariable [namespace current]::time_ini
label $w.top.right.traj.time.steplabel -text "Step:"
entry $w.top.right.traj.time.stepval -width 6 -textvariable [namespace current]::time_step
pack $w.top.right.traj.time.0 $w.top.right.traj.time.inilabel $w.top.right.traj.time.inival $w.top.right.traj.time.steplabel $w.top.right.traj.time.stepval -side left -anchor w
frame $w.top.right.traj.file -relief ridge -bd 0
pack $w.top.right.traj.file -side top -fill x
checkbutton $w.top.right.traj.file.plot -text "Plot" -variable [namespace current]::plot_sw
checkbutton $w.top.right.traj.file.0 -text "Save:" -variable [namespace current]::save_sw\
-command [namespace code {
if {$save_sw} {
$w.top.right.traj.file.name config -state normal
} else {
$w.top.right.traj.file.name config -state disable
}
}]
entry $w.top.right.traj.file.name -width 15 -textvariable [namespace current]::save_file -state disable
pack $w.top.right.traj.file.plot $w.top.right.traj.file.0 $w.top.right.traj.file.name -side left -anchor w
# Iterative Fitting
if {[llength [info procs [namespace current]::iterativeFit]] == 1} {
[namespace current]::iterativeFitGUI
}
# Data
frame $w.data -relief ridge -bd 2
pack $w.data -side top -fill both -expand yes
grid columnconfigure $w.data 1 -weight 1
grid rowconfigure $w.data 1 -weight 1
label $w.data.header_id -text "id" -width 4 -relief sunken
label $w.data.header_mol -text "mol" -width 30 -relief sunken
label $w.data.header_avg -text "avg" -width 7 -relief sunken
label $w.data.header_sd -text "sd" -width 7 -relief sunken
label $w.data.header_min -text "min" -width 7 -relief sunken
label $w.data.header_max -text "max" -width 7 -relief sunken
label $w.data.header_num -text "num" -width 4 -relief sunken
grid $w.data.header_id -column 0 -row 0
grid $w.data.header_mol -column 1 -row 0 -sticky we
grid $w.data.header_avg -column 2 -row 0
grid $w.data.header_sd -column 3 -row 0
grid $w.data.header_min -column 4 -row 0
grid $w.data.header_max -column 5 -row 0
grid $w.data.header_num -column 6 -row 0
set datalist(id) [listbox $w.data.body_id -height 10 -width 4 -relief sunken -exportselection 0 -yscrollcommand [namespace current]::data_yset -selectmode extended]
set datalist(mol) [listbox $w.data.body_mol -height 10 -width 30 -relief sunken -exportselection 0 -yscrollcommand [namespace current]::data_yset -selectmode extended]
set datalist(avg) [listbox $w.data.body_avg -height 10 -width 7 -relief sunken -exportselection 0 -yscrollcommand [namespace current]::data_yset -selectmode extended]
set datalist(sd) [listbox $w.data.body_sd -height 10 -width 7 -relief sunken -exportselection 0 -yscrollcommand [namespace current]::data_yset -selectmode extended]
set datalist(min) [listbox $w.data.body_min -height 10 -width 7 -relief sunken -exportselection 0 -yscrollcommand [namespace current]::data_yset -selectmode extended]
set datalist(max) [listbox $w.data.body_max -height 10 -width 7 -relief sunken -exportselection 0 -yscrollcommand [namespace current]::data_yset -selectmode extended]
set datalist(num) [listbox $w.data.body_num -height 10 -width 4 -relief sunken -exportselection 0 -yscrollcommand [namespace current]::data_yset -selectmode extended]
grid $w.data.body_id -column 0 -row 1 -sticky ns
grid $w.data.body_mol -column 1 -row 1 -sticky nswe
grid $w.data.body_avg -column 2 -row 1 -sticky ns
grid $w.data.body_sd -column 3 -row 1 -sticky ns
grid $w.data.body_min -column 4 -row 1 -sticky ns
grid $w.data.body_max -column 5 -row 1 -sticky ns
grid $w.data.body_num -column 6 -row 1 -sticky ns
foreach key [array names datalist] {
bind $w.data.body_$key <<ListboxSelect>> "[namespace current]::multiple_sel %W"
}
label $w.data.footer_id -text "" -width 4 -anchor e -relief sunken
label $w.data.footer_mol -text "Overall:" -width 30 -anchor e -relief sunken
label $w.data.footer_avg -textvariable [namespace current]::datatot(avg) -width 7 -anchor e -relief sunken
label $w.data.footer_sd -textvariable [namespace current]::datatot(sd) -width 7 -anchor e -relief sunken
label $w.data.footer_min -textvariable [namespace current]::datatot(min) -width 7 -anchor e -relief sunken
label $w.data.footer_max -textvariable [namespace current]::datatot(max) -width 7 -anchor e -relief sunken
label $w.data.footer_num -textvariable [namespace current]::datatot(num) -width 4 -anchor e -relief sunken
grid $w.data.footer_id -column 0 -row 2
grid $w.data.footer_mol -column 1 -row 2 -sticky we
grid $w.data.footer_avg -column 2 -row 2
grid $w.data.footer_sd -column 3 -row 2
grid $w.data.footer_min -column 4 -row 2
grid $w.data.footer_max -column 5 -row 2
grid $w.data.footer_num -column 6 -row 2
# Scrollbar
scrollbar $w.data.scrbar -orient vert -command [namespace current]::data_yview
#scrollbar $w.scrbar.scrbar -relief raised -activerelief raised -bd 2 -elementborderwidth 2 -orient vert -command {rmsdtt::scroll_data}
grid $w.data.scrbar -column 7 -row 0 -rowspan 3 -sticky ns
# Add/remove molecules from the list
frame $w.bottom
pack $w.bottom -side bottom -fill x
button $w.bottom.delall -text "Erase all" -relief raised -command [namespace current]::mol_del
button $w.bottom.del -text "Erase selected" -relief raised -command "[namespace current]::mol_del 1"
button $w.bottom.addall -text "Add all" -relief raised -command [namespace current]::mol_add
button $w.bottom.add -text "Add active" -relief raised -command "[namespace current]::mol_add 1"
pack $w.bottom.delall $w.bottom.del $w.bottom.addall $w.bottom.add -side left -fill x -expand yes
# Final code
[namespace current]::mol_add 1
[namespace current]::update_swap_types
[namespace current]::update_weighted_mol
[namespace current]::ctrlgui
update
wm minsize $w [winfo width $w] [winfo height $w]
wm resizable $w 1 1
}
proc rmsdtt::doRmsd {} {
variable traj_sw
variable traj_all
variable traj_ref
variable save_sw
variable save_file
variable plot_sw
variable colorize
variable rmsd_base
variable swap_sw
variable swap_print
variable rmsd
variable ref_mol
variable ref_frames
variable rms_sel
variable datalist
variable datatot
variable stats
# Parse selection
set rms_sel [[namespace current]::set_sel]
if {$rms_sel == ""} {
showMessage "Selection is empty selection!"
return -code return
}
#puts "DEBUG: rms_sel: $rms_sel"
# Get reference mol/frames
set ref_mol [[namespace current]::get_ref_mol]
if {$rmsd_base == "selected"} {
set sel_index [lindex [$datalist(mol) curselection] 0]
}
if {$rmsd_base == "ave"} {
set ref_frames 0
} else {
set ref_frames {}
if {$traj_sw} {
if {$traj_all} {
for {set i 0} {$i < [molinfo $ref_mol get numframes]} {incr i} {
lappend ref_frames $i
}
} else {
if {$traj_ref >= [molinfo $ref_mol get numframes]} {
showMessage "Frame ref out of range (max is [expr {[molinfo $ref_mol get numframes]-1}])"
return -code return
}
lappend ref_frames $traj_ref
}
} else {
lappend ref_frames [molinfo $ref_mol get frame]
}
}
#puts "DEBUG: ref_mol: $ref_mol"
#puts "DEBUG: ref_frames: $ref_frames"
# Get target mol/frames
# Get all mols to work with
set target_mol [$datalist(id) get 0 end]
#puts "DEBUG: target_mol: $target_mol"
# Calculate average structure
if {$rmsd_base == "ave"} {
set ave_coor [get_ave_coor $target_mol $rms_sel]
#puts $ave_coor
}
# Check number of atoms
if {$rmsd_base == "ave"} {
set ref_natoms [llength $ave_coor]
} else {
set ref_natoms [[atomselect $ref_mol $rms_sel frame 0] num]
}
if {$ref_natoms == 0} {
showMessage "No atoms have been selected!"
return -code return
}
set message ""
foreach i $target_mol {
if {[[atomselect $i $rms_sel frame 0] num] != $ref_natoms } {
append message "$ref_mol ($ref_natoms)\t\t$i ([[atomselect $i $rms_sel frame 0] num])\n"
}
}
if {$message != ""} {
set message "Number of atoms selected differ for molecules:\n$message"
showMessage $message
return -code return
}
# Calculate rmsd and averages
array unset rmsd
array unset rms_ave
array unset rms_val
set rms_tot 0.0
if {$rmsd_base != "ave"} {
set ref_sel [atomselect $ref_mol $rms_sel]
}
foreach i $target_mol {
set target_frames [[namespace current]::get_frames_for_mol $i]
#puts "DEBUG: target_frames($i): $target_frames"
set target_sel [atomselect $i $rms_sel]
set rms_ave($i) 0.0
foreach j $ref_frames {
if {$rmsd_base != "ave"} {
$ref_sel frame $j
}
foreach k $target_frames {
if {$ref_mol == $i && $j == $k} {
continue
}
#puts -nonewline "DEBUG: computing rmsd($ref_mol:$j,$i:$k)"
$target_sel frame $k
if {$rmsd_base == "ave"} {
set rmsd($ref_mol:$j,$i:$k) [get_rmsd_ave $ave_coor $target_sel]
} else {
set rmsd($ref_mol:$j,$i:$k) [get_rmsd $ref_sel $target_sel]
}
#puts " = $rmsd($ref_mol:$j,$i:$k)"
set rms_ave($i) [expr {$rms_ave($i) + $rmsd($ref_mol:$j,$i:$k)}]
}
}
set rms_tot [expr {$rms_tot + $rms_ave($i)}]
set count($i) [llength [array names rmsd *,$i:*]]
#puts "DEBUG: rms_sum($i) = $rms_ave($i) ; count($i) = $count($i)"
if {$count($i) != 0} {
set rms_ave($i) [expr {$rms_ave($i)/$count($i)}]
}
#puts "DEBUG: rms_ave($i) = $rms_ave($i)"
}
set rms_tot [expr {$rms_tot/[array size rmsd]}]
#puts "DEBUG: rms_tot = $rms_tot (count = [array size rmsd])"
# Calculate statistics: standard deviation, min and max
if {$stats} {
set sd_tot 0.0
set random $rmsd([lindex [array names rmsd] 0])
set min_tot $random
set max_tot $random
foreach i $target_mol {
set rms_sd($i) 0.0
if {$count($i) == 0} {
set rms_min($i) 0.0
set rms_max($i) 0.0
continue
}
set random $rmsd([lindex [array names rmsd *,$i:*] 0])
set rms_min($i) $random
set rms_max($i) $random
foreach data [array names rmsd *,$i:*] {
set temp [expr {$rmsd($data) - $rms_ave($i)}]
set rms_sd($i) [expr {$rms_sd($i) + $temp*$temp}]
if {$rmsd($data) < $rms_min($i) } {set rms_min($i) $rmsd($data)}
if {$rmsd($data) > $rms_max($i) } {set rms_max($i) $rmsd($data)}
}
set sd_tot [expr {$sd_tot + $rms_sd($i)}]
if {$count($i) == 1} {
set rms_sd($i) [expr {sqrt($rms_sd($i))}]
} else {
set rms_sd($i) [expr {sqrt( $rms_sd($i) / ($count($i)-1) )}]
}
#puts "DEBUG: rms_sd($i) = $rms_sd($i)"
if {$rms_min($i) < $min_tot} {set min_tot $rms_min($i)}
if {$rms_max($i) > $max_tot} {set max_tot $rms_max($i)}
}
set count_tot [array size rmsd]
if {$count_tot == 1} {
set sd_tot [expr {sqrt($sd_tot)}]
} else {
set sd_tot [expr {sqrt($sd_tot / ($count_tot - 1))}]
}
#puts "DEBUG: sd_tot = $sd_tot"
}
# Reveal values in GUI
foreach v [array names datalist] {
$datalist($v) delete 0 end
}
foreach i $target_mol {
$datalist(id) insert end [format "%d" $i]
$datalist(mol) insert end [format "%s" [molinfo $i get name]]
$datalist(avg) insert end [format "%8.3f" $rms_ave($i)]
if {$stats} {
$datalist(sd) insert end [format "%8.3f" $rms_sd($i)]
$datalist(min) insert end [format "%8.3f" $rms_min($i)]
$datalist(max) insert end [format "%8.3f" $rms_max($i)]
$datalist(num) insert end [format "%4d" $count($i)]
} else {
$datalist(sd) insert end ""
$datalist(min) insert end ""
$datalist(max) insert end ""
$datalist(num) insert end ""
}
}
if {$rmsd_base == "selected"} {
foreach key [array names datalist] {
$datalist($key) selection clear 0 end
$datalist($key) selection set $sel_index
}
}
set datatot(avg) [format "%8.3f" $rms_tot]
if {$stats} {
set datatot(sd) [format "%8.3f" $sd_tot]
set datatot(min) [format "%8.3f" $min_tot]
set datatot(max) [format "%8.3f" $max_tot]
set datatot(num) [format "%4d" [array size rmsd]]
}
if {$traj_sw && $plot_sw && !$traj_all && $colorize} {
[namespace current]::color_data 1
} else {
[namespace current]::color_data 0
}
# Save and plot data
if {$traj_sw && $save_sw} { saveData $save_file }
if {$traj_sw && $plot_sw && !$traj_all} { doPlot }
#puts "DEBUG: ---------------------------------"
# Add selection to history
[namespace current]::ListHistoryPullDownMenu
}
proc rmsdtt::get_ref_mol {} {
variable rmsd_base
variable datalist
switch $rmsd_base {
top {
set mol [molinfo top]
}
ave {
set mol "ave"
}
selected {
set sel_index [lindex [$datalist(mol) curselection] 0]
if {$sel_index == ""} {
showMessage "No molecule has been selected!"
return -code return
}
set mol [$datalist(id) get $sel_index]
}
}
return $mol
}
proc rmsdtt::get_weights { sel } {
variable weighted_mol
variable weighted_field
variable rmsd_base
variable datalist
if {$weighted_mol eq "ref"} {
set mol [[namespace current]::get_ref_mol]
} else {
set mol $weighted_mol
}
return [[atomselect $mol [$sel text]] get $weighted_field]
}
proc rmsdtt::get_rmsd { sel1 sel2 } {
variable swap_sw
variable swap_print
variable weighted_sw
variable weighted_mol
variable weighted_field
if {$weighted_sw} {
set weights [[namespace current]::get_weights $sel1]
set rmsd -1
if [catch { set rmsd [measure rmsd $sel1 $sel2 weight $weights] } msg] {
set message "$msg\nFirst weights were (mol \"$weighted_mol\"; field \"$weighted_field\"):\n[lrange $weights 0 20]..."
showMessage $message
return -code return
}
puts $rmsd
} else {
set rmsd [measure rmsd $sel1 $sel2]
}
if {$swap_sw} {
if {$swap_print} {
puts "\nDEBUG: Swapped residues:"
}
set swapped {}
set mol2 [$sel2 molid]
set frame2 [$sel2 frame]
set sel_text [$sel2 text]
set res [lsort -unique -integer [[atomselect $mol2 "$sel_text and resname [array names ::swap::swap_list]" frame $frame2] get residue]]
foreach r $res {
set s [atomselect $mol2 "residue $r"]
::swap::swap_residue $s $frame2
if {$weighted_sw} {
set rmsd2 [measure rmsd $sel1 $sel2 weight $weights]
} else {
set rmsd2 [measure rmsd $sel1 $sel2]
}
if {$rmsd2 < $rmsd} {
if {$swap_print} {
puts "swapped mol $mol2 frame $frame2 residue $r ([lindex [$s get {resname resid chain segname}] 0]) $rmsd $rmsd2"
}
set rmsd $rmsd2
lappend swapped $s
} else {
::swap::swap_residue $s $frame2
}
}
foreach s $swapped {
if {$s == ""} {
continue
} else {
::swap::swap_residue $s $frame2
}
}
}
return $rmsd
}
proc rmsdtt::get_ave_coor {mols sel_text} {
variable traj_sw
set natoms [[atomselect [lindex $mols 0] $sel_text frame 0] num]
set initialize 1
set ave_coor {}
set tot_frames 0
foreach i $mols {
set nframes [molinfo $i get numframes]
if {$traj_sw && $nframes > 1} {
set avpos [measure avpos [atomselect $i $sel_text] first 0 last [expr {$nframes-1}] step 1]
#puts "$i avpos: $avpos"
if {$initialize} {
for {set j 0} {$j < $natoms} {incr j} {
lappend ave_coor [vecscale [lindex $avpos $j] $nframes]
}
set initialize 0
} else {
for {set j 0} {$j < $natoms} {incr j} {
lset ave_coor $j [vecadd [lindex $ave_coor $j] [vecscale [lindex $avpos $j] $nframes]]
}
}
} else {
set coor [[atomselect $i $sel_text frame [molinfo $i get frame]] get {x y z}]
#puts "$i coor: $coor"
set nframes 1
if {$initialize} {
set ave_coor $coor
set initialize 0
} else {
for {set j 0} {$j < $natoms} {incr j} {
lset ave_coor $j [vecadd [lindex $ave_coor $j] [lindex $coor $j]]
}
}
}
set tot_frames [expr {$tot_frames + $nframes}]
#puts "$i sum: $ave_coor"
}
#puts "tot_frames: $tot_frames"
for {set j 0} {$j < $natoms} {incr j} {
lset ave_coor $j [vecscale [lindex $ave_coor $j] [expr {double(1)/$tot_frames}]]
}
#puts "AVE_COOR: $ave_coor"
return $ave_coor
}
proc rmsdtt::get_rmsd_ave {ref_coor target_sel} {
set target_coor [$target_sel get {x y z}]
set numatoms [llength $ref_coor]
set rmsd 0
# faster than foreach r $ref_coor t $target_sel
for {set k 0} {$k < $numatoms} {incr k} {
set rmsd [expr { $rmsd + [veclength2 [vecsub [lindex $ref_coor $k] [lindex $target_coor $k]]] }]
}
set rmsd [expr {sqrt($rmsd/$numatoms)}]
return $rmsd
}
proc rmsdtt::doAlign {} {
variable w
variable rmsd_base
variable traj_ref
variable traj_sw
variable traj_all
variable datalist
if {$traj_all} {
showMessage "All option not available for Alignment! Deselect it and select a frame reference."
return -code return
}
if {$rmsd_base == "ave"} {
showMessage "Average option not available for Alignment!"
return -code return
}
set ref_mol [[namespace current]::get_ref_mol]
if {$rmsd_base == "selected"} {
set sel_index [lindex [$datalist(mol) curselection] 0]
}
set rms_sel [[namespace current]::set_sel]
set sel_ref [atomselect $ref_mol $rms_sel frame $traj_ref]
set target_mol [$datalist(id) get 0 end]
foreach i $target_mol {
set sel [atomselect $i $rms_sel]
if {$traj_sw == 0} {
if {$i != $ref_mol} {
align $sel [molinfo $i get frame] $sel_ref
}
} else {
for {set j 0} {$j < [molinfo $i get numframes]} {incr j} {
if {$i == $ref_mol && $j == $traj_ref} {
continue
}
align $sel $j $sel_ref
}
}
}
if {$rmsd_base == "selected"} {
foreach key [array names datalist] {
$datalist($key) selection clear 0 end
$datalist($key) selection set $sel_index
}
}
# Add selection to history
[namespace current]::ListHistoryPullDownMenu
}
proc rmsdtt::draw_equiv {} {
variable equiv_sw
variable equiv_byres
variable equiv_atom
variable equiv_lines
variable equiv_mol
variable datalist
[namespace current]::ctrlgui
set equiv_mol [molinfo top]
foreach l $equiv_lines {
if {[graphics $equiv_mol exists $l]} {
graphics $equiv_mol delete $l
}
}
if {!$equiv_sw} {
return
}
set target_mol [$datalist(id) get 0 end]
set sel [set_sel]
lappend equiv_lines [graphics $equiv_mol color 4]
# TODO: too many loops
for {set i 0} {$i < [llength $target_mol]} {incr i} {
set mol1 [lindex $target_mol $i]
if {$equiv_byres} {
set coor1 {}
set residues [lsort -unique -integer [[atomselect $mol1 $sel] get residue]]
foreach r $residues {
set ctmp [[atomselect $mol1 "residue $r and name $equiv_atom"] get {x y z}]
if {[llength $ctmp] == 0} {
set ctmp [[atomselect $mol1 "residue $r"] get {x y z}]
}
lappend coor1 [lindex $ctmp 0]
}
} else {
set coor1 [[atomselect $mol1 $sel] get {x y z}]
}
for {set j [expr {$i+1}]} {$j < [llength $target_mol]} {incr j} {
set mol2 [lindex $target_mol $j]
if {$equiv_byres} {
set coor2 {}
set residues [lsort -unique -integer [[atomselect $mol2 $sel] get residue]]
foreach r $residues {
set ctmp [[atomselect $mol2 "residue $r and name $equiv_atom"] get {x y z}]
if {[llength $ctmp] == 0} {
set ctmp [[atomselect $mol2 "residue $r"] get {x y z}]
}
lappend coor2 [lindex $ctmp 0]
}
} else {
set coor2 [[atomselect $mol2 $sel] get {x y z}]
}
if {[llength $coor1] != [llength $coor2]} {
puts "mol$mol1 and mol$mol2 have different number of atoms"
continue
}
foreach c1 $coor1 c2 $coor2 {
lappend equiv_lines [graphics $equiv_mol line $c1 $c2 width 1 style dashed]
}
}
}
}
proc rmsdtt::align {sel1 frame1 sel2} {
variable weighted_sw
variable weighted_mol
variable weighted_field
$sel1 frame $frame1
if {$weighted_sw} {
set weights [[namespace current]::get_weights $sel1]
set tmatrix ""
if [catch { set tmatrix [measure fit $sel1 $sel2 weight $weights] } msg] {
set message "$msg\nFirst weights were (mol \"$weighted_mol\"; field \"$weighted_field\"):\n[lrange $weights 0 20]..."
showMessage $message
return -code return
}
} else {
set tmatrix [measure fit $sel1 $sel2]
}
set move_sel [atomselect [$sel1 molid] "all" frame $frame1]
$move_sel move $tmatrix
return $tmatrix
}
proc rmsdtt::SaveDataBrowse { {type "data"} } {
variable rmsd
if {![array exists rmsd]} {
showMessage "No data available to save yet!"
return -code return
}
set typeList {
{"Data Files" ".dat .txt .out"}
{"Postscript Files" ".ps"}
{"All files" ".*"}
}
set file [tk_getSaveFile -filetypes $typeList -defaultextension ".dat" -title "Select file to save data" -parent .rmsdtt]
if {$file == ""} {
return
}
if {$type == "data"} {
[namespace current]::saveData $file
} else {
[namespace current]::saveSummary $file
}
}
proc rmsdtt::saveData { file } {
variable traj_sw
variable traj_all
variable time_sw
variable time_ini
variable time_step
variable skip_sw
variable skip_start
variable skip_end
variable skip_steps
variable rmsd
variable ref_mol
variable ref_frames
if {$file == ""} {
showMessage "Filename is missing!"
return -code return
}
if {$skip_sw} {
set ini $skip_start
if {$skip_end == "end"} {
#set end [expr {[llength $list] -1}]
} else {
#set end [lsearch $list $skip_end]
}
set steps [expr {$skip_steps+1}]
} else {
set ini 0
set steps 1
}
# Retrieve mols and frames
set maxframe 0
foreach key [array names rmsd] {
set indices [split $key :,]
set mol [lindex $indices 2]
set frame [lindex $indices 3]
lappend target_mol $mol
lappend target_frames($mol) $frame
if {$frame > $maxframe} {
set maxframe $frame
}
}
set target_mol [lsort -unique -integer $target_mol]
foreach i $target_mol {
set target_frames($i) [lsort -unique -integer $target_frames($i)]
}
if {$traj_sw && $traj_all} {
# Header
if {$time_sw} {
set output "ref_mol\tref_time\tmol\t time\t rmsd\n"
} else {
set output "ref_mol\tref_frame\tmol\tframe\t rmsd\n"
}
foreach k $ref_frames {
set ref_time [expr {$time_ini + $time_step * $k}]
foreach i $target_mol {
foreach j $target_frames($i) {
if {![info exists rmsd($ref_mol:$k,$i:$j)]} {
continue
}
if {$time_sw} {
set time [expr {$time_ini + $time_step * $j}]
append output [format "%7d\t%8.2f\t%3d\t%8.2f\t" $ref_mol $ref_time $i $time]
} else {