# Script reading coordinate files (in any VMD-readable format - NAMDBIN, PDB, DCD, etc.), adding them all as frames in trajectory and saving the trajectory in peaces (to avoid extremely big DCD files)
# Inspired by John Stone's animatepdbs
# Expanded to handle file masks (similar to 'ls' command) without specification of the first and last files in sequence (autodetect). 
# Usage: coordinates2dcd <filesPath> <fileNameFormat> <|fileType> <|outFrames> <|start> <|end>
# For big file sequences can automatically wrap them into trajectory pieces by 'outFrames' frames number.
# Andriy Anishkin (anishkin@icqmail.com) UMCP

proc coordinates2dcd {filesPath fileNameFormat {fileType ""} {outFrames ""} {start ""} {end ""}} {
	if {[llength $filesPath] > 0} {;	# Go to the working folder
		cd $filesPath
	}
	
	# ____ Autodetect the first and last files numbers
	if {([llength $start] == 0)|([llength $end] == 0)} {;	#start or end numbers were not specified
		set fileNameMask [regsub {%[0-9]*[dui]{1}} $fileNameFormat "*"]
		set fileList [lsort [ls $fileNameMask]]
		set start [scan [lindex $fileList 1] $fileNameFormat];	# not the 0 element - this is the parent folder
		set end [scan [lindex $fileList end] $fileNameFormat]
		puts "[expr {[llength $fileList]-1}] file(s) found: from [lindex $fileList 1]   to   [lindex $fileList end]"
	}
	
	# ____ Load the first file
	set fileName [format $fileNameFormat [expr {$start}]]
	puts "Reading initial file $fileName in the sequence"
	if {[llength $fileType]>0} {
		mol new $fileName type $fileType waitfor all
	} else {
		mol new $fileName waitfor all
	}
	
	# ____ Load the rest of files
	set firstDcdFrame 0
	set filesRead 1
	puts "Reading coordinates files and writing sets of frames into DCD..."
	incr start
	for {set i $start} {$i <= $end} {incr i 1} {
		set fileName [format $fileNameFormat [expr {$i}]]
		
		if {[llength $fileType]>0} {
			mol addfile $fileName type $fileType waitfor all
		} else {
			mol addfile $fileName waitfor all
		}
		
		incr filesRead
		 
		set framesRead [molinfo top get numframes]
		if {[llength $outFrames]>0} {
			if {$framesRead>=$outFrames} {;	# Maximal number of frames in the output file reached
				set lastDcdFrame [expr {$firstDcdFrame+$framesRead-1}]
				set firstDcdFrameString [format "%0[string length $end]d" $firstDcdFrame]
				set lastDcdFrameString [format "%0[string length $end]d" $lastDcdFrame]
				animate write dcd [join "[molinfo top get name] _frames $firstDcdFrameString to $lastDcdFrameString .dcd" {}] beg 0 end [expr {$outFrames-1}] waitfor all 
				set firstDcdFrame [expr {$lastDcdFrame+1}]
				animate delete beg 0 end [expr {$outFrames-1}]
				
			}
		}
	}
	
	
	if {[molinfo top get numframes]>0} {
		set lastDcdFrame [expr {$firstDcdFrame+[molinfo top get numframes]-1}]
		set firstDcdFrameString [format "%0[string length $end]d" $firstDcdFrame]
		set lastDcdFrameString [format "%0[string length $end]d" $lastDcdFrame]
		animate write dcd [join "[molinfo top get name] _frames $firstDcdFrameString to $lastDcdFrameString .dcd" {}] waitfor all 
		animate delete all
	}
	
	puts "[expr {$lastDcdFrame+1}] frames from $filesRead files read"
	
	puts "Writing DCD trajectory finished!"
}