diff --git a/.project b/.project new file mode 100644 index 0000000..e626399 --- /dev/null +++ b/.project @@ -0,0 +1,11 @@ + + + QC_JSON_Schema + + + + + + + + diff --git a/Coding/Java/src/org/jmol/adapters/readers/quantum/QCJSONReader.java b/Coding/Java/src/org/jmol/adapters/readers/quantum/QCJSONReader.java new file mode 100644 index 0000000..adf3a32 --- /dev/null +++ b/Coding/Java/src/org/jmol/adapters/readers/quantum/QCJSONReader.java @@ -0,0 +1,377 @@ +package org.jmol.adapter.readers.quantum; + +import java.util.ArrayList; +import java.util.Arrays; +import java.util.Hashtable; +import java.util.Map; + +import javajs.util.AU; +import javajs.util.Lst; +import javajs.util.SB; + +import org.jmol.adapter.smarter.Atom; +import org.jmol.api.JmolAdapter; +import org.jmol.util.Logger; +import org.qcschema.QCSchemaUnits; + +/** + * A molecular structure and orbital reader for MolDen files. + * See http://www.cmbi.ru.nl/molden/molden_format.html + * + * updated by Bob Hanson for Jmol 12.0/12.1 + * + * adding [spacegroup] [operators] [cell] [cellaxes] for Jmol 14.3.7 + * + * @author Matthew Zwier + */ + +public class QCJSONReader extends MoldenReader { + + private Map job; + + private int jobCount; + + private int modelCount; + + @SuppressWarnings("unchecked") + @Override + protected void initializeReader() { + super.initializeReader(); + SB sb = new SB(); + try { + while (rd() != null) + sb.append(line); + Lst json = vwr.parseJSONArray(sb.toString()); + // first record is version tag + Logger.info(json.get(0).toString()); + // second record is Jmol info; not used here + jobCount = json.size() - 2; + for (int i = 0; i < jobCount; i++) + processJob((Map)json.get(i + 2)); + } catch (Exception e) { + e.printStackTrace(); + } + continuing = false; + } + + /** + * @param job + * @throws Exception + */ + private void processJob(Map job) throws Exception { + this.job = job; + readSteps(); + /* + if (loadVibrations) + readFreqsAndModes(); + if (loadGeometries) + readGeometryOptimization(); + checkSymmetry(); + if (asc.atomSetCount == 1 && moData != null) + finalizeMOData(moData); + */ + } + + @Override + public void finalizeSubclassReader() throws Exception { + finalizeReaderASCR(); + } + + private void readSteps() throws Exception { + ArrayList steps = QCSchemaUnits.getList(job, "steps"); + int nSteps = steps.size(); + for (int iStep = 0; iStep < nSteps; iStep++) { + if (!doGetModel(++modelCount, null)) { + if (!checkLastModel()) + return; + continue; + } + asc.newAtomSet(); + @SuppressWarnings("unchecked") + Map step = (Map) steps.get(iStep); + Map topology = getMapSafely(step, "topology"); + Map atoms = getMapSafely(topology, "atoms"); + + // one or the other of these is required: + String[] symbols = QCSchemaUnits.getStringArray(atoms, "symbol"); + int[] atomNumbers = QCSchemaUnits.getIntArray(atoms, "atom_number"); + String[] atom_names = QCSchemaUnits.getStringArray(atoms, "atom_names"); + + double[] coords = QCSchemaUnits.getDoubleArray(atoms, "coords"); + modelAtomCount = coords.length / 3; + double f = QCSchemaUnits.getConversionFactor(atoms, "coords", QCSchemaUnits.UNITS_ANGSTROMS); + boolean isFractional = (f == 0); + setFractionalCoordinates(isFractional); + if (isFractional) { + f = QCSchemaUnits.getConversionFactor(atoms, "unit_cell", QCSchemaUnits.UNITS_ANGSTROMS); + double[] cell = QCSchemaUnits.getDoubleArray(atoms, "unit_cell"); + // a b c alpha beta gamma + // m.m00, m.m10, m.m20, // Va + // m.m01, m.m11, m.m21, // Vb + // m.m02, m.m12, m.m22, // Vc + // dimension, (float) volume, + if (cell == null) { + Logger.error("topology.unit_cell is missing even though atoms are listed as fractional"); + } else { + for (int i = 0; i < 6; i++) { + switch (i) { + case 3: + f = 1; + //$FALL-THROUGH$ + default: + setUnitCellItem(i, (float)(cell[i] * f)); + break; + } + } + } + } + for (int i = 0, pt = 0; i < modelAtomCount; i++) { + Atom atom = asc.addNewAtom(); + setAtomCoordXYZ(atom, (float)(coords[pt++] * f), (float)(coords[pt++] * f), (float) (coords[pt++] + * f)); + String sym = (symbols == null ? JmolAdapter + .getElementSymbol(atomNumbers[i]) : symbols[i]); + atom.atomName = (atom_names == null ? sym : atom_names[i]); + atom.elementNumber = (short) (atomNumbers == null ? JmolAdapter + .getElementNumber(sym) : atomNumbers[i]); + } + if (doReadMolecularOrbitals) { + readMolecularOrbitals(getMapSafely(step, "molecular_orbitals")); + clearOrbitals(); + } + applySymmetryAndSetTrajectory(); + if (loadVibrations) { + readFreqsAndModes(QCSchemaUnits.getList(step, "vibrations")); + } + + } + } + + private boolean readFreqsAndModes(ArrayList vibrations) throws Exception { + // "frequency":{"value":-0.00,"units":["cm^-1","?"]}, + // "ir_intensity":{"value":0.000005,"units":["au",1]}, + // "vectors":[ + + if (vibrations != null) { + int n = vibrations.size(); + for (int i = 0; i < n; i++) { + @SuppressWarnings("unchecked") + Map vib = (Map) vibrations.get(i); + double freq = QCSchemaUnits.getDouble(vib, "frequency", QCSchemaUnits.UNITS_CM_1); + double[] vectors = QCSchemaUnits.getDoubleArray(vib, "vectors"); + if (i > 0) + asc.cloneLastAtomSet(); + asc.setAtomSetFrequency(null, null, "" + freq, QCSchemaUnits.UNITS_CM_1); + int i0 = asc.getLastAtomSetAtomIndex(); + for (int j = 0, pt = 0; j < modelAtomCount; j++) { + asc.addVibrationVector(j + i0, (float) (vectors[pt++] * ANGSTROMS_PER_BOHR), + (float) (vectors[pt++] * ANGSTROMS_PER_BOHR), (float) (vectors[pt++] + * ANGSTROMS_PER_BOHR)); + } + } + } + return true; + } + + private boolean haveEnergy = true; + + /** + * Read basis and orbital information. + * + * @param molecular_orbitals + * @return true if successful + * + * @throws Exception + */ + private boolean readMolecularOrbitals(Map molecular_orbitals) throws Exception { + if (molecular_orbitals == null) + return false; + String moBasisID = molecular_orbitals.get("basis_id").toString();//:"MOBASIS_1" + if (!readBasis(moBasisID)) + return false; + Boolean isNormalized = (Boolean) molecular_orbitals.get("__jmol_normalized"); + if (isNormalized != null && isNormalized.booleanValue()) + moData.put("isNormalized", isNormalized); + calculationType = (String) molecular_orbitals.get("__jmol_calculation_type"); + if (calculationType == null) + calculationType = "?"; + moData.put("calculationType", calculationType); + + ArrayList mos = QCSchemaUnits.getList(molecular_orbitals, "orbitals"); + int n = mos.size(); + for (int i = 0; i < n; i++) { + @SuppressWarnings("unchecked") + Map thisMO = (Map) mos.get(i); + double energy = QCSchemaUnits.getDouble(thisMO, "energy", "ev"); + double occupancy = QCSchemaUnits.getDouble(thisMO, "occupancy", null); + String symmetry = (String) thisMO.get("symmetry"); + String spin = (String) thisMO.get("type"); + if (spin != null) { + if (spin.indexOf("beta") >= 0) + alphaBeta = "beta"; + else if (spin.indexOf("alpha") >= 0) + alphaBeta = "alpha"; + } + float[] coefs = toFloatArray(QCSchemaUnits.getDoubleArray(thisMO, "coefficients")); + line = "" + symmetry; + if (filterMO()) { + Map mo = new Hashtable(); + mo.put("coefficients", coefs); + if (Double.isNaN(energy)) { + haveEnergy = false; + } else { + mo.put("energy", Float.valueOf((float) energy)); + } + if (!Double.isNaN(occupancy)) + mo.put("occupancy", Float.valueOf((float) occupancy)); + if (symmetry != null) + mo.put("symmetry", symmetry); + if (alphaBeta.length() > 0) + mo.put("type", alphaBeta); + setMO(mo); + if (debugging) { + Logger.debug(coefs.length + " coefficients in MO " + orbitals.size()); + } + } + } + if (debugging) + Logger.debug("read " + orbitals.size() + " MOs"); + ArrayList units = QCSchemaUnits.getList(molecular_orbitals, "orbitals_energy_units"); + String sunits = (units == null ? null : units.get(0).toString()); + setMOs(sunits == null || sunits.equals("?") ? "?" : sunits); + if (haveEnergy && doSort) + sortMOs(); + return false; + } + + private float[] toFloatArray(double[] da) { + float[] fa = new float[da.length]; + for (int j = da.length; --j >= 0;) + fa[j] = (float) da[j]; + return fa; + } + + String lastBasisID = null; + private boolean readBasis(String moBasisID) throws Exception { + Map moBasisData = getMapSafely(job, "mo_bases"); + Map moBasis = getMapSafely(moBasisData, moBasisID); + if (moBasis == null) { + Logger.error("No job.mo_bases entry for " + moBasisID); + return false; + } + if (moBasisID == lastBasisID) + return true; + lastBasisID = moBasisID; + ArrayList listG = QCSchemaUnits.getList(moBasis, "gaussians"); + ArrayList listS = QCSchemaUnits.getList(moBasis, "shells"); + if (listG == null && listS == null) { + listG = listS = QCSchemaUnits.getList(moBasis, "slaters"); + } + if ((listG == null) != (listS == null)) { + Logger.error("gaussians/shells or slaters missing"); + return false; + } + if (listG == listS) { + readSlaterBasis(listS); + } else { + readGaussianBasis(listG, listS); + } + return true; + } + + boolean readSlaterBasis(ArrayList listS) throws Exception { + /* + 1 0 0 0 1 1.5521451600 0.9776767193 + 1 1 0 0 0 1.5521451600 1.6933857512 + 1 0 1 0 0 1.5521451600 1.6933857512 + 1 0 0 1 0 1.5521451600 1.6933857512 + 2 0 0 0 0 1.4738648100 1.0095121222 + 3 0 0 0 0 1.4738648100 1.0095121222 + */ + + nCoef = listS.size(); + for (int i = 0; i < nCoef; i++) { + double[] a = QCSchemaUnits.getDoubleArray(listS.get(i), null); + addSlater((int) a[0], (int) a[1], (int) a[2], (int) a[3], (int) a[4], (float) a[5], (float) a[6]); + } + setSlaters(false, false); + return true; + } + + private boolean readGaussianBasis(ArrayList listG, ArrayList listS) throws Exception { + shells = new Lst(); + for (int i = 0; i < listS.size(); i++) + shells.addLast(QCSchemaUnits.getIntArray(listS.get(i), null)); + int gaussianPtr = listG.size(); + float[][] garray = AU.newFloat2(gaussianPtr); + // [[exp, coef], [exp, coef],...] with sp [exp, coef1, coef2] + for (int i = 0; i < gaussianPtr; i++) + garray[i] = toFloatArray(QCSchemaUnits.getDoubleArray(listG.get(i), null)); + moData.put("shells", shells); + moData.put("gaussians", garray); + Logger.info(shells.size() + " slater shells read"); + Logger.info(garray.length + " gaussian primitives read"); + //Logger.info(nCoef + " MO coefficients expected for orbital type " + orbitalType); + asc.setCurrentModelInfo("moData", moData); + return false; + } + + @SuppressWarnings("unchecked") + private void sortMOs() { + Object[] list = orbitals.toArray(new Object[orbitals.size()]); + Arrays.sort(list, new MOEnergySorter()); + orbitals.clear(); + for (int i = 0; i < list.length; i++) + orbitals.addLast((Map)list[i]); + } + + /** + * Safely get a Map from a Map using a key. + * @param map + * @param key + * @return the Map or null + */ + @SuppressWarnings("unchecked") + private static Map getMapSafely(Map map, String key) { + return (map == null ? null : (Map) map.get(key)); + } + + /////////////////// from Molden reader -- TODO ///////////////// + +// private boolean checkSymmetry() throws Exception { +// // extension for symmetry +// if (line.startsWith("[SPACEGROUP]")) { +// setSpaceGroupName(rd()); +// rd(); +// return true; +// } +// if (line.startsWith("[OPERATORS]")) { +// while (rd() != null && line.indexOf("[") < 0) +// if (line.length() > 0) { +// Logger.info("adding operator " + line); +// setSymmetryOperator(line); +// } +// return true; +// } +// if (line.startsWith("[CELL]")) { +// rd(); +// Logger.info("setting cell dimensions " + line); +// // ANGS assumed here +// next[0] = 0; +// for (int i = 0; i < 6; i++) +// setUnitCellItem(i, parseFloat()); +// rd(); +// return true; +// } +// if (line.startsWith("[CELLAXES]")) { +// float[] f = new float[9]; +// fillFloatArray(null, 0, f); +// addExplicitLatticeVector(0, f, 0); +// addExplicitLatticeVector(1, f, 3); +// addExplicitLatticeVector(2, f, 6); +// return true; +// } +// return false; +// } +// +} diff --git a/Coding/Java/src/org/jmol/adapters/writers/QCJSONWriter.java b/Coding/Java/src/org/jmol/adapters/writers/QCJSONWriter.java new file mode 100644 index 0000000..408810f --- /dev/null +++ b/Coding/Java/src/org/jmol/adapters/writers/QCJSONWriter.java @@ -0,0 +1,609 @@ +package org.jmol.adapter.writers; + +import java.io.OutputStream; +import java.util.Date; +import java.util.Hashtable; +import java.util.Map; + +import javajs.util.DF; +import javajs.util.Lst; +import javajs.util.P3; +import javajs.util.PT; +import javajs.util.SB; + +import org.jmol.api.SymmetryInterface; +import org.jmol.java.BS; +import org.jmol.modelset.Atom; +import org.jmol.quantum.SlaterData; +import org.jmol.util.JSONWriter; +import org.jmol.util.Vibration; +import org.jmol.viewer.Viewer; +import org.qcschema.QCSchemaUnits; + +/** + * A very experimental class for writing QCJSON files. This standard is in the + * process of being developed, so any of this could change at any time. + * + * All we have here is Bob Hanson's experiment with getting Jmol to save and + * restore structures, vibrations, and molecular orbitals. + * + * Data set Bob is using is at + * + * https://sourceforge.net/p/jmol/code/HEAD/tree/trunk/Jmol-datafiles/qcjson + * + */ +public class QCJSONWriter extends JSONWriter { + + // Current status: + // + // 2017.12.14 Generating valid JSON code that can be read back in for tested files. + // + + private Map moBases = new Hashtable(); + + private boolean filterMOs; + + private Viewer vwr; + + public void set(Viewer viewer, OutputStream os) { + vwr = viewer; + setWriteNullAsString(false); + setStream(os); + } + + @Override + public String toString() { + return (oc == null ? "{}" : oc.toString()); + } + + public void writeJSON() { + openSchema(); + writeMagic(); + oc.append(",\n"); + writeSchemaMetadata(); + writeJobs(); + closeSchema(); + } + + public void writeSchemaMetadata() { + mapOpen(); + mapAddKeyValue("__jmol_created", new Date(), ",\n"); + mapAddKeyValue("__jmol_source", vwr.getP("_modelFile"), ""); + mapClose(); + } + + public void openSchema() { + arrayOpen(false); + } + + public void writeMagic() { + writeString(QCSchemaUnits.version); + } + + public void closeSchema() { + oc.append("\n"); + arrayClose(false); + closeStream(); + } + + public void writeJobs() { + // only one job in Jmol + writeJob(1); + } + + public void writeJob(int iJob) { + append(",\n"); + mapOpen(); + { + mapAddKeyValue("__jmol_block", "Job " + iJob, ",\n"); + writeJobMetadata(); + writeModels(); + writeMOBases(); + } + mapClose(); + } + + public void writeJobMetadata() { + mapAddKey("metadata"); + mapOpen(); + { + mapAddMapAllExcept("__jmol_info", vwr.getModelSetAuxiliaryInfo(), + ";group3Counts;properties;group3Lists;models;unitCellParams;"); + } + mapClose(); + } + + public void writeModels() { + int nModels = vwr.ms.mc; + oc.append(",\n"); + mapAddKey("steps"); + arrayOpen(true); + { + oc.append("\n"); + for (int i = 0; i < nModels;) { + if (i > 0) + append(",\n"); + i = writeModel(i); + } + } + arrayClose(true); + } + + public int writeModel(int modelIndex) { + int nextModel = modelIndex + 1; + append(""); + mapOpen(); + { + mapAddKeyValue("__jmol_block", "Model " + (modelIndex + 1), ",\n"); + writeTopology(modelIndex); + if (isVibration(modelIndex)) { + oc.append(",\n"); + nextModel = writeVibrations(modelIndex); + } + if (haveMOData(modelIndex)) { + oc.append(",\n"); + writeMOData(modelIndex); + } + oc.append(",\n"); + writeModelMetadata(modelIndex); + } + mapClose(); + oc.append("\n"); + return nextModel; + } + + public void writeTopology(int modelIndex) { + mapAddKey("topology"); + mapOpen(); + { + writeAtoms(modelIndex); + writeBonds(modelIndex); + } + mapClose(); + } + + public Object getProperty(int modelIndex, String key) { + @SuppressWarnings("unchecked") + Map props = (Map) (modelIndex >= vwr.ms.am.length ? null + : vwr.ms.am[modelIndex].auxiliaryInfo.get("modelProperties")); + return (props == null ? null : props.get(key)); + } + + private boolean isVibration(int modelIndex) { + return (vwr.ms.getLastVibrationVector(modelIndex, 0) >= 0); + } + + public void writeModelMetadata(int modelIndex) { + mapAddKey("metadata"); + mapOpen(); + { + mapAddMapAllExcept("__jmol_info", vwr.ms.am[modelIndex].auxiliaryInfo, + ";.PATH;PATH;fileName;moData;unitCellParams;"); + } + mapClose(); + } + + public void writeAtoms(int modelIndex) { + SparseArray symbols = new SparseArray("_RLE_"); + SparseArray numbers = new SparseArray("_RLE_"); + SparseArray charges = new SparseArray("_RLE_"); + SparseArray names = new SparseArray("_RLE_"); + SparseArray types = new SparseArray("_RLE_"); + mapAddKey("atoms"); + mapOpen(); + { + SymmetryInterface unitCell = vwr.ms.getUnitCell(modelIndex); + boolean isFractional = (unitCell != null && !unitCell.isBio()); + if (isFractional) { + float[] params = unitCell.getUnitCellAsArray(false); + writePrefix_Units("unit_cell", "angstroms"); + mapAddKeyValue("unit_cell", params, ",\n"); + } + writePrefix_Units("coords", isFractional ? "fractional" : "angstroms"); + mapAddKey("coords"); + arrayOpen(true); + { + oc.append("\n"); + BS bs = vwr.getModelUndeletedAtomsBitSet(modelIndex); + int last = bs.length() - 1; + P3 pt = new P3(); + for (int i = bs.nextSetBit(0); i >= 0; i = bs.nextSetBit(i + 1)) { + Atom a = vwr.ms.at[i]; + append(""); + pt.setT(a); + if (isFractional) + unitCell.toFractional(pt, false); + oc.append(formatNumber(pt.x)).append(",\t") + .append(formatNumber(pt.y)).append(",\t") + .append(formatNumber(pt.z)).append(i < last ? ",\n" : "\n"); + symbols.add(PT.esc(a.getElementSymbol())); + numbers.add("" + a.getElementNumber()); + charges.add("" + a.getPartialCharge()); + String name = a.getAtomName(); + names.add(name); + String type = a.getAtomType(); + types.add(type.equals(name) ? null : type); + } + } + arrayClose(true); + oc.append(",\n"); + if (charges.isNumericAndNonZero()) { + mapAddKeyValueRaw("charge", charges, ",\n"); + } + if (types.hasValues()) { + mapAddKeyValueRaw("types", types, ",\n"); + } + mapAddKeyValueRaw("symbol", symbols, ",\n"); + mapAddKeyValueRaw("atom_number", numbers, "\n"); + } + mapClose(); + } + + private String formatNumber(float x) { + return (x < 0 ? "" : " ") + DF.formatDecimal(x, -6); + } + + private void writePrefix_Units(String prefix, String units) { + mapAddKeyValueRaw(prefix + "_units", QCSchemaUnits.getUnitsJSON(units, false), + ",\n"); + } + + public void writeBonds(int modelIndex) { + // TODO + } + + public int writeVibrations(int modelIndex) { + mapAddKey("vibrations"); + arrayOpen(true); + { + oc.append("\n"); + String sep = null; + int ivib = 0; + modelIndex--; + while (isVibration(++modelIndex)) { + if (sep != null) + oc.append(sep); + sep = ",\n"; + append(""); + mapOpen(); + { + mapAddKeyValue("__jmol_block", "Vibration " + (++ivib), ",\n"); + Object value = getProperty(modelIndex, "FreqValue"); + String freq = (String) getProperty(modelIndex, "Frequency"); + String intensity = (String) getProperty(modelIndex, "IRIntensity"); + String[] tokens; + if (value == null) { + System.out.println("model " + modelIndex + + " has no _M.properties.FreqValue"); + } + if (freq == null) { + System.out.println("model " + modelIndex + + " has no _M.properties.Frequency"); + } else { + tokens = PT.split(freq, " "); + if (tokens.length == 1) { + System.out.println("model " + modelIndex + + " has no frequency units"); + } + writeMapKeyValueUnits("frequency", value, tokens[1]); + } + if (intensity != null) { + tokens = PT.split(intensity, " "); + writeMapKeyValueUnits("ir_intensity", tokens[0], tokens[1]); + + } + String label = (String) getProperty(modelIndex, "FrequencyLabel"); + if (label != null) + mapAddKeyValue("label", label, ",\n"); + mapAddKey("vectors"); + arrayOpen(true); + { + oc.append("\n"); + BS bs = vwr.getModelUndeletedAtomsBitSet(modelIndex); + int last = bs.length() - 1; + for (int i = bs.nextSetBit(0); i >= 0; i = bs.nextSetBit(i + 1)) { + Atom a = vwr.ms.at[i]; + Vibration v = a.getVibrationVector(); + append(""); + oc.append(formatNumber(v.x)).append(",\t") + .append(formatNumber(v.y)).append(",\t") + .append(formatNumber(v.z)).append(i < last ? ",\n" : "\n"); + } + } + arrayClose(true); + } + append(""); + mapClose(); + } + } + oc.append("\n"); + arrayClose(true); + return modelIndex; + } + + private void writeMapKeyValueUnits(String key, Object value, String units) { + mapAddKeyValueRaw(key, "{\"value\":" + value + ",\"units\":" + + QCSchemaUnits.getUnitsJSON(units, false) + "}", ",\n"); + } + + private boolean haveMOData(int modelIndex) { + return (getAuxiliaryData(modelIndex, "moData") != null); + } + + private Object getAuxiliaryData(int modelIndex, String key) { + return vwr.ms.am[modelIndex].auxiliaryInfo.get(key); + } + + private int basisID = 0; + private Lst shells; + + private int[][] dfCoefMaps; + + private void writeMOData(int modelIndex) { + @SuppressWarnings("unchecked") + Map moData = (Map) getAuxiliaryData( + modelIndex, "moData"); + Map moDataJSON = new Hashtable(); + moDataJSON.put("orbitals", moData.get("mos")); + // units + String units = (String) moData.get("EnergyUnits"); + if (units == null) + units = "?"; + moDataJSON.put("orbitals_energy_units", QCSchemaUnits.getUnitsJSON(units, true)); + // normalization is critical for Molden, NWChem, and many other readers. + // not needed for Gaussian, Jaguar, WebMO, Spartan, or GenNBO + moDataJSON.put("__jmol_normalized", + Boolean.valueOf(moData.get("isNormalized") == Boolean.TRUE)); + String type = (String) moData.get("calculationType"); + moDataJSON.put("__jmol_calculation_type", type == null ? "?" : type); + // @SuppressWarnings("unchecked") + // Map orbitalMaps = (Map) moData.get("orbitalMaps"); + // if (orbitalMaps != null && !orbitalMaps.isEmpty()) { + // moDataJSON.put("jmol_orbital_maps", orbitalMaps); + // } + moDataJSON.put("basis_id", getBasisID(moData)); + filterMOs = true; + setModifyKeys(fixIntegration()); + mapAddKeyValue("molecular_orbitals", moDataJSON, "\n"); + setModifyKeys(null); + filterMOs = false; + append(""); + } + + private static Map integrationKeyMap; + + /** + * When an MO is calculated in Jmol, Jmol will check the integration so that + * it can be checked to be close to 1.0000. This integration value is saved + * back in the MO data, but it is not a standard key. (As though anything is + * here!) + * + * So we set a key mapping to replace it. + * + * @return the "integration" key map + */ + private static Map fixIntegration() { + if (integrationKeyMap == null) { + integrationKeyMap = new Hashtable(); + integrationKeyMap.put("integration", "__jmol_integration"); + } + return integrationKeyMap; + } + + @Override + protected Object getAndCheckValue(Map map, String key) { + if (filterMOs) { + if (key.equals("dfCoefMaps")) + return null; + if (key.equals("symmetry")) + return ((String) map.get(key)).replace('_', ' ').trim(); + if (key.equals("coefficients") && dfCoefMaps != null) { + return fixCoefficients((double[]) map.get(key)); + } + + } + return map.get(key); + } + + /** + * Jmol allows for a set of arrays that map coefficient indicies with + * nonstandard order to Gaussian/Molden order. Here we do the conversion upon + * writing so that the order is always Gaussian/Molden order. + * + * @param coeffs + * @return + */ + private Object fixCoefficients(double[] coeffs) { + double[] c = new double[coeffs.length]; + for (int i = 0, n = shells.size(); i < n; i++) { + int[] shell = shells.get(i); + int type = shell[1]; + int[] map = dfCoefMaps[type]; + for (int j = 0, coefPtr = 0; j < map.length; j++, coefPtr++) + c[coefPtr + j] = coeffs[coefPtr + map[j]]; + } + return c; + } + + @SuppressWarnings("unchecked") + private String getBasisID(Map moData) { + String hash = "!"; + dfCoefMaps = (int[][]) moData.get("dfCoefMaps"); + if (dfCoefMaps != null) { + // just looking for a non-zero map + boolean haveMap = false; + for (int i = 0; !haveMap && i < dfCoefMaps.length; i++) { + int[] m = dfCoefMaps[i]; + for (int j = 0; j < m.length; j++) + if (m[j] != 0) { + haveMap = true; + break; + } + } + if (!haveMap) + dfCoefMaps = null; + } + Object gaussians = moData.get("gaussians"); + if (gaussians != null) { + hash += gaussians.hashCode(); + } + shells = (Lst) moData.get("shells"); + if (shells != null) { + hash += shells.hashCode(); + } + Object slaters = moData.get("slaters"); + if (slaters != null) { + hash += slaters.hashCode(); + } + String key = (String) moBases.get(hash); + if (key == null) { + moBases.put(hash, key = "MOBASIS_" + ++basisID); + Map map = new Hashtable(); + if (gaussians != null) { + map.put("gaussians", gaussians); + } + if (shells != null) { + + // shells array: [iAtom, type, gaussianPtr, gaussianCount] + // + // where type is one of: + // + // final public static int S = 0; + // final public static int P = 1; + // final public static int SP = 2; + // final public static int DS = 3; + // final public static int DC = 4; + // final public static int FS = 5; + // final public static int FC = 6; + // final public static int GS = 7; + // final public static int GC = 8; + // final public static int HS = 9; + // final public static int HC = 10; + // final public static int IS = 11; + // final public static int IC = 12; + + // Note that this is currently implemented in Jmol with reference to a + // coefficient map that allows us to maintain the file-based MO ordering + // and only map the actual coefficient to the function at MO creation time. + + map.put("shells", shells); + } + if (slaters != null) { + map.put("slaters", slaters); + } + moBases.put(key, map); + } + return key; + } + + public void writeMOBases() { + if (moBases.isEmpty()) + return; + oc.append(",\n"); + mapAddKey("mo_bases"); + mapOpen(); + { + String sep = ""; + for (String key : moBases.keySet()) { + if (key.startsWith("!")) + continue; + append(sep); + mapAddKeyValue(key, moBases.get(key), "\n"); + sep = ","; + } + } + mapClose(); + moBases.clear(); + } + + @Override + public void writeObject(Object o) { + if (o instanceof SlaterData) { + oc.append(o.toString()); + } else { + super.writeObject(o); + } + } + + //// sparse array handling //// + public class SparseArray extends SB { + private int repeatCount = 0; + private int elementCount = 0; + private String lastElement = null; + private String sep = ""; + private String type; // _RLE_ + private boolean isRLE; + + public SparseArray(String type) { + this.type = type; + isRLE = (type.equals("_RLE_")); + } + + protected void add(String element) { + if (element == null) + element = "null"; + if (!isRLE) { + append(sep); + append(element); + sep = ","; + return; + } + if (repeatCount > 0 && !element.equals(lastElement)) { + append(sep); + appendI(repeatCount); + sep = ","; + append(sep); + append(lastElement); + repeatCount = 0; + } + lastElement = element; + repeatCount++; + elementCount++; + } + + public String lastElement() { + return lastElement; + } + + public boolean isEmpty() { + return (elementCount == 0); + } + + public boolean allNaN() { + return (allSame() && PT.parseFloat(lastElement) == Float.NaN); + } + + public boolean allNull() { + return (allSame() && lastElement.equals("null")); + } + + public boolean allEmptyString() { + return (allSame() && lastElement.equals("")); + } + + public boolean allSame() { + return (!isEmpty() && elementCount == repeatCount); + } + + public boolean allZero() { + return (allSame() && PT.parseFloat(lastElement) != Float.NaN); + } + + public boolean hasValues() { + return (!allSame() || !allNull() && !allEmptyString()); + } + + public boolean isNumericAndNonZero() { + return (allSame() && !allNaN() && !allZero()); + } + + @Override + public String toString() { + String s = super.toString(); + return (s.length() == 0 ? "[]" : "[\"" + type + "\"," + s + + (repeatCount > 0 ? sep + repeatCount + "," + lastElement : "") + + "]"); + } + } + +} diff --git a/Coding/Java/src/org/qcschema/QCSchemaUnits.java b/Coding/Java/src/org/qcschema/QCSchemaUnits.java new file mode 100644 index 0000000..0d5846f --- /dev/null +++ b/Coding/Java/src/org/qcschema/QCSchemaUnits.java @@ -0,0 +1,358 @@ +package org.qcschema; + +import java.util.Hashtable; +import java.util.Map; + +import java.util.ArrayList; + +import org.jmol.viewer.Viewer; + +/** + * A general Java class for working with QCShema units and array types. + * + * j2sNative blocks can be ignored -- they just increase efficiency in the JavaScript rendition of Jmol. + * + */ +public class QCSchemaUnits { + + public final static String version = "QCJSON 0-0-0.Jmol_" + + Viewer.getJmolVersion().replace(' ', '_'); + + // + // source: http://cccbdb.nist.gov/hartree.asp + // A hartree is equal to 2625.5 kJ/mol, 627.5 kcal/mol, 27.211 eV, and 219474.6 cm-1. + // One bohr = 0.529 177 210 67 x 10-10 m + + public final static String UNITS_FRACTIONAL = "fractional"; + + public final static String UNITS_AU = "au"; + public final static double TOAU_AU = 1; + + // distance + + public final static String UNITS_CM = "cm"; + public final static double TOAU_CM = 1/0.52917721067e-8; + + public final static String UNITS_M = "m"; + public final static double TOAU_M = 1/0.52917721067e-10; + + public final static String UNITS_ANGSTROMS = "angstroms"; + public final static double TOAU_ANGSTROMS = 1/0.52917721067; // 1.88972613; + + public final static String UNITS_BOHR = "bohr"; + public final static double TOAU_BOHR = 1; + + // energy + + public final static String UNITS_HARTREE = "hartree"; + public final static double TOAU_HARTREE = 1; + + public final static String UNITS_EV = "ev"; + public final static double TOAU_EV = 1/27.211; //0.03688675765; + + public final static String UNITS_CM_1 = "cm-1"; + public final static double TOAU_CM_1 = 1/219474.6; //4.5563359e-6; + + public final static String UNITS_KJ_MOL = "kj/mol"; + public final static double TOAU_KJ_MOL = 1/2635.5; //0.00038087983; + + public final static String UNITS_KCAL_MOL = "kcal/mol"; + public final static double TOAU_KCAL_MOL = 1/627.5; //0.00159362549; + + /** + * A very simple and efficient way to catalog string matches. Far faster than ENUM. + * Note that singular or plural on anstroms, bohrs, or hartrees both work. + */ + private final static String knownUnits = + /////0 1 2 3 4 5 6 7 8 + /////012345678901234567890123456789012345678901234567890123456789012345678901234567890123 + "cm cm^-1 cm-1 angstroms au atomic units fractional bohrs hartrees ev kj_mol kcal_mol"; + + private static Hashtable htConvert = new Hashtable(); + + /** + * Get the standard conversion factor to atomic units for this unit. + * + * @param units + * @return the nominal conversion factor or 0 ("fractional") or Double.NaN (unknown) + */ + public static double getFactorToAU(String units) { + switch (knownUnits.indexOf(units.toLowerCase())) { + case 0: + // units = UNITS_CM + return TOAU_CM; + case 1: + // units = UNITS_M + return TOAU_M; + case 3: + case 9: + //units = UNITS_CM_1; + return TOAU_CM_1; + case 14: + //units = UNITS_ANGSTROMS; + return TOAU_ANGSTROMS; + case 24: + case 27: + //units = UNITS_AU; + return 1; + case 40: + //units = "UNITS_FRACTIONAL"; + return 0; + case 51: + //units = UNITS_BOHR; + return TOAU_BOHR; + case 57: + //units = UNITS_HARTREE; + return TOAU_HARTREE; + case 66: + //units = UNITS_EV; + return TOAU_EV; + case 69: + //units = UNITS_KCAL_MOL; + return TOAU_KCAL_MOL; + case 76: + //units = UNITS_KJ_MOL; + return TOAU_KJ_MOL; + default: + return Double.NaN; + } + } + + /** + * Calculate the unit conversion between two units, using a static + * unit-to-unit cache for efficiency. + * + * Not used in Jmol. + * + * @param fromUnits + * @param toUnits + * @return conversion factor or Double.NaN if anything goes wrong. + */ + public static double getUnitConversion(String fromUnits, String toUnits) { + if (fromUnits.equalsIgnoreCase(toUnits)) + return 1; + String key = "" + fromUnits + toUnits; + Double d = htConvert.get(key); + if (d != null) + return d.doubleValue(); + double val = Double.NaN; + try { + double toAUDesired = getFactorToAU(toUnits); + double toAUActual = getFactorToAU(fromUnits); + val = toAUActual / toAUDesired; + } catch (Exception e) { + // just leave it as 1 + } + htConvert.put(key, Double.valueOf(val)); + return val; + } + + /** + * For a reader, use the JSON [units, factor] along with a desired unit to get the conversion + * factor from file values to desired units. + * + * Currently, this method only looks at the factor in the JSON if we do not already know the conversion factor. + * + * @param unitsFactor [units, factor] list or null if to AU is desired. + * @param unitsDesired + * @return the conversion factor or Double.NaN if not uncodable + */ + public static double getConversionFactorTo(ArrayList unitsFactor, String unitsDesired) { + try { + double toAUDesired = getFactorToAU(unitsDesired); + double toAUActual = getFactorToAU(unitsFactor == null ? UNITS_AU : unitsFactor.get(0).toString()); + if (Double.isNaN(toAUActual)) + toAUActual = Double.parseDouble(unitsFactor.get(1).toString()); + return toAUActual / toAUDesired; + } catch (Exception e) { + return Double.NaN; + } + } + + /** + * Read a {value:xxxx, units:["name",toAU]} map, converting it to the desired units. + * + * @param valueUnits + * @param toUnits + * @return converted value + */ + public static double convertValue(Map valueUnits, String toUnits) { + return getDouble(valueUnits, "value", null) * getConversionFactor(valueUnits, "units", toUnits); + } + + /** + * Get the [name, toAU] JSON code or just a new String[] {name, toAU}. + * If the conversion is not known, return [name, "?"] + * @param name + * @param asArray + * @return String or String[] + */ + public static Object getUnitsJSON(String name, boolean asArray) { + double d = getFactorToAU(name); + String toAU = (!Double.isNaN(d) ? "" + d : asArray ? "?" : "\"?\""); + return (asArray ? new String[] { name, toAU } : "[\"" + name + "\"," + + toAU + "]"); + } + + + /** + * Get the necessary conversion factor to the desired units from a key_units or atomic units + * @param map + * @param key map key that has associated key_units element or null for "from atomic units" + * @param toUnits + * @return conversion factor + */ + public static double getConversionFactor(Map map, String key, String toUnits) { + ArrayList list = getList(map, key + "_units"); + String units = (list == null ? null : list.get(0).toString()); + double f = getConversionFactorTo(list, toUnits); + if (Double.isNaN(f)) { + System.out.println("units for " + units + "? " + units); + f = 1; + } + return f; + } + + /** + * Reads a value from an associative array, converting it to the desired units. + * + * @param map + * @param key + * @param toUnits + * @return value + */ + @SuppressWarnings("unchecked") + public static double getDouble(Map map, String key, String toUnits) { + Object o = map.get(key); + double conv = 1; + if (toUnits != null) + if (o instanceof Map) { + // "frequency":{"value":-0.00,"units":["cm^-1",4.5563359e-6]}, + return convertValue((Map) o, toUnits); + } else if (map.containsKey(key + "_units")) { + // "frequency_units":["cm^-1",4.5563359e-6], + // "frequency":-0.00, + conv = getConversionFactor(map, key, toUnits); + } + return (o == null ? Double.NaN : ((Number) o).doubleValue() * conv); + } + + /** + * Retrieve an array of any sort as a list of objects, possibly unpacking it + * if it is run-length encoded. + * + * @param mapOrList + * @param key + * @return unpacked array + */ + public static ArrayList getList(Object mapOrList, String key) { + @SuppressWarnings("unchecked") + ArrayList list = (ArrayList) (key == null ? mapOrList + : ((Map) mapOrList).get(key)); + if (list == null) + return null; + int n = list.size(); + if (n == 0 || !"_RLE_".equals(list.get(0))) + return list; + ArrayList list1 = newList(); + for (int i = 1; i < n; i++) { + int count = ((Number) list.get(i)).intValue(); + Object value = list.get(++i); + for (int j = 0; j < count; j++) + /** + * j2s avoids overloaded add() for speed + * + * @j2sNative list1.addLast(value); + */ + { + list1.add(value); + } + } + return list1; + } + + /** + * @return ArrayList, or in JavaScript javajs.util.Lst + * + * @j2sNative + * + * return new javajs.util.Lst(); + */ + protected static ArrayList newList() { + return new ArrayList(); + } + + /** + * Retrieve a double array, possibly unpacking it if it is run-length encoded. + * Read any error as Double.NaN. + * + * @param mapOrList + * @param key into mapOrList, or null if mapOrList is a list + * @return unpacked double[] + */ + public static double[] getDoubleArray(Object mapOrList, String key) { + ArrayList list = getList(mapOrList, key); + if (list == null) + return null; + double[] a = new double[list.size()]; + for (int i = a.length; --i >= 0;) { + try { + a[i] = ((Number) list.get(i)).doubleValue(); + } catch (Exception e) { + a[i] = Double.NaN; + } + } + return a; + } + + /** + * Retrieve an int array, possibly unpacking it if it is run-length encoded. + * Any error causes this method to return null. + * + * @param mapOrList the list to unpack, or map to pull the list form using the key + * @param key the map key, or null if mapOrList is already a list + * @return unpacked int[] or null if mapOrList is null or there is an error + */ + public static int[] getIntArray(Object mapOrList, String key) { + ArrayList list = getList(mapOrList, key); + if (list != null) { + try { + int[] a = new int[list.size()]; + for (int i = a.length; --i >= 0;) + a[i] = ((Number) list.get(i)).intValue(); + return a; + } catch (Exception e) { + // return null in this case + } + } + return null; + } + + /** + * Retrieve a String array, possibly unpacking it if it is run-length encoded. + * Any "null" string is read as null. + * + * @param mapOrList the list to unpack, or map to pull the list form using the key + * @param key the map key, or null if mapOrList is already a list + * @return unpacked string[] or null if mapOrList is null + */ + public static String[] getStringArray(Object mapOrList, String key) { + ArrayList list = getList(mapOrList, key); + if (list == null) + return null; + String[] a = new String[list.size()]; + for (int i = a.length; --i >= 0;) { + Object o = list.get(i); + a[i] = (o == null ? null : list.get(i).toString()); + } + return a; + } + +// static { +// System.out.println(getUnitConversion("Angstroms", "cm")); +// System.out.println(getUnitConversion("Angstroms", "bohr")); +// System.out.println(getUnitConversion("bohr", "Angstroms")); +// System.out.println(getUnitConversion("AX", "bohr")); +// System.out.println(getUnitConversion("bohr", "AX")); +// } +} diff --git a/Examples/28ce209867afd272d361a00322960160 b/Examples/28ce209867afd272d361a00322960160 deleted file mode 160000 index b4bb403..0000000 --- a/Examples/28ce209867afd272d361a00322960160 +++ /dev/null @@ -1 +0,0 @@ -Subproject commit b4bb4030240e71c2c37879fcd945efef254c58f5 diff --git a/Examples/NWChemOutputToJson b/Examples/NWChemOutputToJson deleted file mode 160000 index 427d799..0000000 --- a/Examples/NWChemOutputToJson +++ /dev/null @@ -1 +0,0 @@ -Subproject commit 427d7991fc91221ddd0e39cd60845b2252cb982f diff --git a/Examples/chemicaljson b/Examples/chemicaljson deleted file mode 160000 index 416c10d..0000000 --- a/Examples/chemicaljson +++ /dev/null @@ -1 +0,0 @@ -Subproject commit 416c10dcb97aaee2d5dfa26b7403f117944e8287 diff --git a/Examples/json_schema b/Examples/json_schema deleted file mode 160000 index cd72ef3..0000000 --- a/Examples/json_schema +++ /dev/null @@ -1 +0,0 @@ -Subproject commit cd72ef3e156090d9e16127ed3aa4ac2295c7dc44 diff --git a/Examples/molecular-design-toolkit.wiki b/Examples/molecular-design-toolkit.wiki deleted file mode 160000 index d70bbfe..0000000 --- a/Examples/molecular-design-toolkit.wiki +++ /dev/null @@ -1 +0,0 @@ -Subproject commit d70bbfebb5d731516fd08f2af26bd1ef2f18d2ac