From: toshinagata1964 <toshinagata1964@a2be9bc6-48de-4e38-9406-05402d4bc13c>
Date: Mon, 26 Sep 2011 15:50:15 +0000 (+0000)
Subject: Improve crd/mdcrd import.
X-Git-Tag: v1.0.2~511
X-Git-Url: http://git.osdn.net/view?a=commitdiff_plain;h=bbda4277d98a3f4de5d4ccf90e29358fd0d8acf4;p=molby%2FMolby.git

Improve crd/mdcrd import.

git-svn-id: svn+ssh://svn.sourceforge.jp/svnroot/molby/trunk@122 a2be9bc6-48de-4e38-9406-05402d4bc13c
---

diff --git a/MolLib/MD/MDCore.c b/MolLib/MD/MDCore.c
index 7bf2eef..64ea1a4 100644
--- a/MolLib/MD/MDCore.c
+++ b/MolLib/MD/MDCore.c
@@ -1406,9 +1406,16 @@ md_prepare(MDArena *arena, int check_only)
 				err = "cannot create log file";
 		}
 		if (err == NULL && arena->coord_result_name != NULL && arena->coord_result == NULL) {
+			char molname[64];
 			arena->coord_result = fopen(arena->coord_result_name, "wb");
 			if (arena->coord_result == NULL)
 				err = "cannot create coord file";
+			else {
+				int natoms = (arena->output_expanded_atoms ? arena->mol->natoms : arena->natoms_uniq);
+				MoleculeCallback_displayName(mol, molname, sizeof molname);
+				fprintf(arena->coord_result, "TITLE: %s, %d atoms\n", molname, natoms);
+				fflush(arena->coord_result);
+			}
 		}
 		if (err == NULL && arena->vel_result_name != NULL && arena->vel_result == NULL) {
 			arena->vel_result = fopen(arena->vel_result_name, "wb");
@@ -1914,13 +1921,16 @@ md_output_results(MDArena *arena)
 				VecScaleInc(r, *bv, ap[i].wrap_dy);
 				VecScaleInc(r, *cv, ap[i].wrap_dz);
 			}
-			fprintf(arena->coord_result, " %7.3f%s %7.3f%s %7.3f%s",
+			fprintf(arena->coord_result, "%8.3f%s%8.3f%s%8.3f%s",
 				r.x, (j == 9 ? "\n" : ""),
 				r.y, (j == 8 ? "\n" : ""),
 				r.z, (j == 7 ? "\n" : ""));
 		}
 		if (j != 0)
 			fprintf(arena->coord_result, "\n");
+		if (arena->periodic_a || arena->periodic_b || arena->periodic_c) {
+			fprintf(arena->coord_result, "%8.3f%8.3f%8.3f\n", arena->mol->cell->cell[0], arena->mol->cell->cell[1], arena->mol->cell->cell[2]);
+		}
 		fflush(arena->coord_result);
 	}
 	if (arena->vel_result != NULL) {
diff --git a/Scripts/loadsave.rb b/Scripts/loadsave.rb
index 0600c15..7828349 100755
--- a/Scripts/loadsave.rb
+++ b/Scripts/loadsave.rb
@@ -25,22 +25,24 @@ class Molecule
     count = 0
 	frame = 0
 	coords = (0...natoms).collect { Vector3D[0, 0, 0] }
-    periodic = (self.box && self.box[4].all? { |n| n != 0 })
+    periodic = (self.box && self.box[4].any? { |n| n != 0 })
 	show_progress_panel("Loading AMBER crd file...")
 #    puts "sframe = #{sframe}, pos = #{fp.pos}"
-    line = fp.gets   #  Skip first line
-    while 1
-      line = fp.gets
-      if line == nil
-	    if (count > 0)
-		  raise MolbyError, sprintf("crd format error - file ended in the middle of frame %d", frame)
-		end
-	    fp.close
-		break
-      end
+    while line = fp.gets
       line.chomp!
-#      values = line.split(' ')
 	  values = line.scan(/......../)
+	  if fp.lineno == 1
+	    #  The first line should be skipped. However, if this line seems to contain
+		#  coordinates, then try reading them
+		if values.size != 10 && (values.size > 10 || values.size != natoms * 3)
+		  next  #  Wrong number of coordinates
+		end
+		if values.each { |v| Float(v) rescue break } == nil
+		  #  The line contains non-number 
+		  next
+	    end
+		puts "Loadcrd: The title line seems to be missing"
+	  end
       if count + values.size > natoms * 3
         raise MolbyError, sprintf("crd format error - too many values at line %d in file %s; number of atoms = %d, current frame = %d", fp.lineno, fp.path, natoms, frame)
       end
@@ -59,8 +61,11 @@ class Molecule
 		if periodic
 		  #  Should have box information
 		  line = fp.gets
-		  if line == nil || (values = line.chomp.split(' ')).length != 3
-		    raise "The molecule has a periodic cell but the crd file does not contain cell information"
+		  if line == nil || (values = line.chomp.scan(/......../)).length != 3
+		    #  Periodic but no cell information
+			puts "Loadcrd: the molecule has a periodic cell but the crd file does not contain cell info"
+			periodic = false
+			redo
 	      end
 		  self.cell = [Float(values[0]), Float(values[1]), Float(values[2]), 90, 90, 90]
 		end
@@ -71,6 +76,10 @@ class Molecule
 		end
       end
     end
+	fp.close
+	if (count > 0)
+	  raise MolbyError, sprintf("crd format error - file ended in the middle of frame %d", frame)
+	end
 	hide_progress_panel
 	if frame > 0
 	  self.frame = self.nframes - 1