//
/*
 * io_ambertraj.cpp
 *
 *  Created on: Apr 18, 2019
 *      Author: Shailesh
 */

#include "io.h"
#include <vector>

#define DIM 3
#define GROMACS_MAGIC 1993
#define NM2ANG 10.

using namespace std;

delphi_integer CIO::readMDCRD(ifstream& ifMdcrdFileStream,
		SPrmtopHeadPointers prmtopHeadPointer, vector<CFrame>& vctOfFrames,
		vector<delphi_real>& vctAtomRadii, const int& iFrameFirst,
		const int& iFrameLast, const int& iFrameStride) {
	int valuesPerFrame = (prmtopHeadPointer.natoms * 3);
	int linesPerFrame = valuesPerFrame / 10;
	if (prmtopHeadPointer.natoms * 3 % 10 != 0)
		linesPerFrame += 1;
	// if box information is present it spans 1-line
	// if (prmtopHeadPointer.ifbox)
	//	linesPerFrame += 1;

	int firstFrameOffset = 1; // trajectory title spans first line of mdcrd file
	int currentLineIndex = 0;

	string strLine;
	getline(ifMdcrdFileStream, strLine);
	int currentFrameIndex = 0;
	while (strLine) {
		bool bRecordFrame = shouldRecordFrame(currentFrameIndex, iFrameFirst,
				iFrameLast, iFrameStride);
		if (bRecordFrame) {
			CFrame newFrame(prmtopHeadPointer.natoms);
			int valueIndex = 0;
			int valuesToRead = valuesPerFrame;
			for (int i = 0; i < linesPerFrame; ++i) {
				getline(ifMdcrdFileStream, strLine);
				int valuesInThisLine = (valuesToRead > 10) ? 10 : valuesToRead;
				for (int j = 0; j < valuesInThisLine; ++j, ++valueIndex) {
					int atomIdx = valueIndex / 3;
					int axisIdx = valueIndex % 3;
					double value = atof(
							removeSpace(strLine.substr(j * 8, 8)).c_str());
					switch (axisIdx) {
					case 0:
						newFrame.coordinates[atomIdx].nX = value;
						break;
					case 1:
						newFrame.coordinates[atomIdx].nY = value;
						break;
					case 2:
						newFrame.coordinates[atomIdx].nZ = value;
						break;
					}
					valuesToRead--;
				}
			}
			newFrame.computeFrameGeometry(vctAtomRadii);
			//newFrame.dumpFrame();
			vctOfFrames.push_back(newFrame);
		} else {
			// skip this frame
			for (int i = 0; i < linesPerFrame; ++i) {
				getline(ifMdcrdFileStream, strLine);
			}
		}
		// if box information is present skipt it
		if (prmtopHeadPointer.ifbox) {
			getline(ifMdcrdFileStream, strLine);
		}
#ifdef VERBOSE
		if (bRecordFrame)
		cout << "io_trr> Reading frame : " << iFrame << endl;
#endif
		currentFrameIndex++;
		if (currentFrameIndex == iFrameLast) {
			break;
		}
	}
	return prmtopHeadPointer.natoms;
}
