1#!/usr/local/bin/perl -w 2 3# in: input filename 4$in = shift || die("Please specify input filename"); 5# If your exchange was every N ps and you saved every M ps you can make for 6# the missing frames by setting extra to (N/M - 1). If N/M is not integer, 7# you're out of luck and you will not be able to demux your trajectories at all. 8$extra = shift || 0; 9$ndx = "replica_index.xvg"; 10$temp = "replica_temp.xvg"; 11 12@comm = ("-----------------------------------------------------------------", 13 "Going to read a file containing the exchange information from", 14 "your mdrun log file ($in).", 15 "This will produce a file ($ndx) suitable for", 16 "demultiplexing your trajectories using trjcat,", 17 "as well as a replica temperature file ($temp).", 18 "Each entry in the log file will be copied $extra times.", 19 "-----------------------------------------------------------------"); 20for($c=0; ($c<=$#comm); $c++) { 21 printf("$comm[$c]\n"); 22} 23 24# Open input and output files 25open (IN_FILE,"$in") || die ("Cannot open input file $in"); 26open (NDX,">$ndx") || die("Opening $ndx for writing"); 27open (TEMP,">$temp") || die("Opening $temp for writing"); 28 29 30sub pr_order { 31 my $t = shift; 32 my $nrepl = shift; 33 printf(NDX "%-20.2f",$t); 34 for(my $k=0; ($k<$nrepl); $k++) { 35 my $oo = shift; 36 printf(NDX " %3d",$oo); 37 } 38 printf(NDX "\n"); 39} 40 41sub pr_revorder { 42 my $t = shift; 43 my $nrepl = shift; 44 printf(TEMP "%-20.2f",$t); 45 for(my $k=0; ($k<$nrepl); $k++) { 46 my $oo = shift; 47 printf(TEMP " %3d",$oo); 48 } 49 printf(TEMP "\n"); 50} 51 52$nrepl = 0; 53$init = 0; 54$tstep = 0; 55$nline = 0; 56$tinit = 0; 57while ($line = <IN_FILE>) { 58 chomp($line); 59 60 if (index($line,"init_t") >= 0) { 61 @log_line = split (' ',$line); 62 $tinit = $log_line[2]; 63 } 64 if (index($line,"Repl") == 0) { 65 @log_line = split (' ',$line); 66 if (index($line,"There") >= 0) { 67 $nrepl = $log_line[3]; 68 } 69 elsif (index($line,"time") >= 0) { 70 $tstep = $log_line[6]; 71 } 72 elsif ((index($line,"Repl ex") == 0) && ($nrepl == 0)) { 73 # Determine number of replicas from the exchange information 74 printf("%s\n%s\n", 75 "WARNING: I did not find a statement about number of replicas", 76 "I will try to determine it from the exchange information."); 77 for($k=2; ($k<=$#log_line); $k++) { 78 if ($log_line[$k] ne "x") { 79 $nrepl++; 80 } 81 } 82 } 83 if (($init == 0) && ($nrepl > 0)) { 84 printf("There are $nrepl replicas.\n"); 85 86 @order = (); 87 @revorder = (); 88 for($k=0; ($k<$nrepl); $k++) { 89 $order[$k] = $k; 90 $revorder[$k] = $k; 91 } 92 for($ee=0; ($ee<=$extra); $ee++) { 93 pr_order($tinit+$ee,$nrepl,@order); 94 pr_revorder($tinit+$ee,$nrepl,@revorder); 95 $nline++; 96 } 97 $init = 1; 98 } 99 100 if (index($line,"Repl ex") == 0) { 101 $k = 0; 102 for($m=3; ($m<$#log_line); $m++) { 103 if ($log_line[$m] eq "x") { 104 $revorder[$order[$k]] = $k+1; 105 $revorder[$order[$k+1]] = $k; 106 $tmp = $order[$k]; 107 $order[$k] = $order[$k+1]; 108 $order[$k+1] = $tmp; 109# printf ("Swapping %d and %d on line %d\n",$k,$k+1,$line_number); 110 } 111 else { 112 $k++; 113 } 114 } 115 for($ee=0; ($ee<=$extra); $ee++) { 116 pr_order($tstep+$ee,$nrepl,@order); 117 pr_revorder($tstep+$ee,$nrepl,@revorder); 118 $nline++; 119 } 120 } 121 } 122} 123close IN_FILE; 124close NDX; 125close TEMP; 126 127printf ("Finished writing $ndx and $temp with %d lines\n",$nline); 128