package name.mjw.jamber.IO.AMBER;
import name.mjw.jamber.IO.FortranFormat;
import org.apache.log4j.Logger;
import javax.vecmath.Point3d;
import javax.vecmath.Vector3d;
import java.io.BufferedReader;
import java.io.FileReader;
import java.io.IOException;
import java.text.ParseException;
import java.util.ArrayList;
import java.util.Iterator;
import java.util.List;
/**
* AMBER INPut CooRDinate
* <p>
*
* This is a representation of the AMBER coordinate file <a
* href="http://ambermd.org/formats.html#restart">format</a>.
*
* @author mjw
*/
public class Inpcrd {
private final Logger LOG = Logger.getLogger(Inpcrd.class);
public Inpcrd(String fileName) {
try {
read(fileName);
} catch (ParseException | IOException e) {
// TODO Auto-generated catch block
e.printStackTrace();
}
}
public Inpcrd() {
}
/**
* The title of the current run, from the AMBER parameter/topology file
*/
private String title;
/**
* Total number of atoms in coordinate file
*/
private int numberOfAtoms = 0;
/**
* Current time in the simulation in picoseconds
*/
private double timeStep = 0.0;
/**
* Atom coordinates in Angstrom
*/
private ArrayList<Point3d> coordinates = new ArrayList<>();
/**
* Atom velocities
*/
private ArrayList<Vector3d> velocities = new ArrayList<>();
/**
* The inpcrd file can contain velocities, set if this true
*/
private Boolean hasVelocities = false;
/**
* The inpcrd file can contain box information, set if this true
*/
private Boolean hasBoxInformation = false;
public String getTitle() {
return title;
}
public void setTitle(String title) {
this.title = title;
}
public int getNumberOfAtoms() {
return numberOfAtoms;
}
public void setNumberOfAtoms(int numberOfAtoms) {
this.numberOfAtoms = numberOfAtoms;
}
public double getTimeStep() {
return timeStep;
}
public void setTimeStep(double timeStep) {
this.timeStep = timeStep;
}
public ArrayList<Point3d> getPositions() {
return coordinates;
}
public void setPositions(ArrayList<Point3d> coordinates) {
this.coordinates = coordinates;
}
public ArrayList<Vector3d> getVelocities() {
return velocities;
}
public void setVelocities(ArrayList<Vector3d> velocities) {
this.velocities = velocities;
}
public Boolean getHasVelocities() {
return hasVelocities;
}
void setHasVelocities(Boolean hasVelocities) {
this.hasVelocities = hasVelocities;
}
public Boolean getHasBoxInformation() {
return hasBoxInformation;
}
void setHasBoxInformation(Boolean hasBoxInformation) {
this.hasBoxInformation = hasBoxInformation;
}
/**
* Given a filename, reads in a Inpcrd file
*
* @param fileName
* Filename of inpcrd file
* @throws ParseException
* @throws IOException
*/
public void read(String fileName) throws ParseException, IOException {
String line;
ArrayList<Object> lineObjects;
try {
BufferedReader in = new BufferedReader(new FileReader(fileName));
line = in.readLine();
// Title
lineObjects = FortranFormat.read(line, "(20A4)");
LOG.debug("read(): " + line);
for (Object lineObject : lineObjects) {
if (lineObject != null) {
this.title = ((String) lineObject);
}
}
// Number of Atoms, time
// 22 0.5005000E+00
line = in.readLine();
lineObjects = FortranFormat.read(line, "(I7,5E15.7)");
LOG.debug("read(): " + line);
this.numberOfAtoms = (Integer) lineObjects.get(0);
// We could have a crd file, which has time information as well
if (lineObjects.get(1) != null) {
this.timeStep = (Double) lineObjects.get(1);
} else {
this.timeStep = 0.0;
}
// Atom coordinates
// 2.0773821 1.0381474 0.1406479 2.1129473 2.1410015 0.1174182
// Read everything in first and then process it
List<Double> rawEntries = new ArrayList<>();
while ((line = in.readLine()) != null) {
lineObjects = FortranFormat.read(line, "(6F12.7)");
LOG.debug("read(): " + line);
for (Object object : lineObjects) {
rawEntries.add((Double) object);
}
}
in.close();
// Work out what we have in this file
LOG.debug("read(): rawEntries.size() is " + rawEntries.size());
if (rawEntries.size() == numberOfAtoms * 3) {
LOG.debug("read(): Only Atom coordinates found");
}
if (rawEntries.size() == ((numberOfAtoms * 3) + 6)) {
LOG.debug("read(): Atom coordinates and boxInformation found");
setHasBoxInformation(true);
}
if (rawEntries.size() == ((numberOfAtoms * 3) * 2)) {
LOG.debug("read(): Atom coordinates and velocities found");
setHasVelocities(true);
}
if (rawEntries.size() == (((numberOfAtoms * 3) * 2) + 6)) {
LOG.debug("read(): Atom coordinates, velocities and boxInformationfound");
setHasBoxInformation(true);
setHasVelocities(true);
}
// TODO what happens if none of the above applies?
// Parse out atom coordinates
List<Double> atomEntries = rawEntries.subList(0,
(numberOfAtoms * 3));
for (Iterator<Double> iterator = atomEntries.iterator(); iterator
.hasNext();) {
// Get x y and z
Double x = iterator.next();
Double y = iterator.next();
Double z = iterator.next();
Point3d coordinate = new Point3d(x, y, z);
LOG.debug(coordinate);
this.coordinates.add(coordinate);
}
if (hasVelocities) {
// Parse out velocities
List<Double> velocityEntries = rawEntries.subList(
(numberOfAtoms * 3), (numberOfAtoms * 3) * 2);
for (Iterator<Double> iterator = velocityEntries.iterator(); iterator
.hasNext();) {
// Get x y and z
Double vx = iterator.next();
Double vy = iterator.next();
Double vz = iterator.next();
Vector3d velocity = new Vector3d(vx, vy, vz);
LOG.debug(velocity);
this.velocities.add(velocity);
}
}
if (hasBoxInformation) {
int end = rawEntries.size();
int start = end - 6;
LOG.debug("Looking for box coordinates starting at: " + start);
LOG.debug("Looking for box coordinates ending at: " + end);
// TODO
// List<Double> boxCoordinates = rawEntries.subList(start, end);
}
} catch (ParseException | IOException e) {
// TODO Auto-generated catch block
e.printStackTrace();
}
}
}