diff --git a/common/python/nomadcore/md_data_access/MDDataAccess.py b/common/python/nomadcore/md_data_access/MDDataAccess.py
index 773bbd124064a717d2f42a154dd3994dad371330..715854a262b927492bb9e847a376b46736a71f47 100644
--- a/common/python/nomadcore/md_data_access/MDDataAccess.py
+++ b/common/python/nomadcore/md_data_access/MDDataAccess.py
@@ -320,8 +320,8 @@ def pymConvertTopoDict(topoStor, topoPYM):
if noninter:
noninter = False
typeList.extend(list([
- atmlist[0][2],
- atmlist[1][2]
+ atom_type_dict[atmlist[0][1]],
+ atom_type_dict[atmlist[1][1]]
]))
interTypeDict.update({interNum : typeList})
interNum += 1
@@ -674,8 +674,8 @@ def mdaConvertTopoDict(topoStor, topoMDA):
if noninter:
noninter = False
typeList.extend(list([
- atmlist[0][2],
- atmlist[1][2]
+ atom_type_dict[atmlist[0][1]],
+ atom_type_dict[atmlist[1][1]]
]))
interTypeDict.update({interNum : typeList})
interNum += 1
@@ -2269,3 +2269,194 @@ if __name__ == "__main__":
atom = MDdata.iread()
print(i)
+ if "mdanalysis" in MDdata.topocode:
+
+ def getseg(seg):
+ regex = re.compile(r"^seg_[0-9]+_b")
+ if regex.findall(seg.segid):
+ segname = bytes(re.sub(r"^seg_[0-9]+_b","",seg.segid).replace("'",""), "utf-8").decode("utf-8")
+ else:
+ segname = seg.segid
+ return [seg.ix,segname]
+
+ def get_atomseg(atom):
+ return [atom.segment.ix,atom.segment.segname.decode('utf-8')]
+
+ def getres(atom):
+ return [atom.resname.decode('utf-8'),atom.resid]
+
+ def getatom(atom):
+ atmid = atom.ix
+ atmname = atom.name.decode('utf-8')
+ atmtyp = atom.type.decode('utf-8')
+ atmres = atom.resname.decode('utf-8')
+ atmresid = atom.resid
+ atmsegid = atom.segid
+ return [atmid,atmname,atmtyp,atmres,atmresid,atmsegid]
+
+ def residuepairs(atom1,atom2):
+ return [getres(atom1),getres(atom2)]
+
+ def atompairs(atom1,atom2):
+ return [getatom(atom1),getatom(atom2)]
+
+ print(MDdata.topohandler.atoms)
+ attrlist = dir(MDdata.topohandler.atoms)
+ if "names" in attrlist:
+ print("Atom Names:",MDdata.topohandler.atoms.names)
+ if "types" in attrlist:
+ print("Atom Types:",MDdata.topohandler.atoms.types)
+ if "masses" in attrlist:
+ print("Atom Masses:",MDdata.topohandler.atoms.masses)
+ if "charges" in attrlist:
+ print("Atom Charges:",MDdata.topohandler.atoms.charges)
+ if "radiuses" in attrlist:
+ print("Atom Radiuses:",MDdata.topohandler.atoms.radiuses)
+ if "bfactors" in attrlist:
+ print("Atom bFactors:",MDdata.topohandler.atoms.bfactors)
+ print("Atom List:",[mda.topology.guessers.guess_atom_element(getatom(a)[1]) for a in MDdata.topohandler.atoms])
+ print("Atom Residues:",MDdata.topohandler.atoms.residues)
+ print("Atom Residue ids:",MDdata.topohandler.atoms.residues.ix)
+ print("Atom Residue Names:",[a.resname for a in MDdata.topohandler.atoms.residues])
+ print("Atom Segments:",MDdata.topohandler.atoms.segments)
+ print("Atom Bonds:",MDdata.topohandler.atoms.bonds)
+ print("Atom Angles:",MDdata.topohandler.atoms.angles)
+ print("Atom Dihedrals:",MDdata.topohandler.atoms.dihedrals)
+ print("Atom Impropers:",MDdata.topohandler.atoms.impropers)
+ print(MDdata.topohandler.bonds)
+ topd = MDdata.topohandler.bonds.topDict
+ topdk = topd.keys()
+ for key in topdk:
+ topt = topd[key]
+ for b in topt:
+ reslist=residuepairs(b.atoms[0],b.atoms[1])
+ atmlist=atompairs(b.atoms[0],b.atoms[1])
+ if str(reslist[0][0]) != str(reslist[1][0]):
+ print("RES:",reslist,
+ "ATM:",atmlist)
+
+
+ #print("Bonds:",[ [a[0].ix,a[1].ix] for a in MDdata.topohandler.bonds])
+ print(MDdata.topohandler.angles)
+ #print("Angles:",[ [a[0].ix,a[1].ix,a[2].ix] for a in MDdata.topohandler.angles])
+ print(MDdata.topohandler.residues)
+ #print("Residue bonds:",[a.bonds for a in MDdata.topohandler.residues])
+ #for r in MDdata.topohandler.residues:
+ # print(r.resname,r.resid)
+ print([getseg(seg)[1] for seg in MDdata.topohandler.segments])
+ #print("Segment Names:",[a.segment for a in MDdata.topohandler.segments])
+ #print("Residue Names:",MDdata.topohandler.residues.names)
+ #print("Residue Types:",MDdata.topohandler.residues.names)
+ #print("Residue Bonds:",MDdata.topohandler.residues.bonds)
+ #print("Residue Angles:",MDdata.topohandler.residues.angles)
+
+
+ print(mda.topology.guessers.guess_atom_element("BERK"))
+
+ #angles = mdtraj.compute_angles(mytopology)
+ #dihedrals = mdtraj.compute_dihedrals(mytopology)
+ #print(dihedrals)
+ #print(' Chains: %s' % [res for res in mytopology.chains])
+ #print(' Residues: %s' % [res for res in mytopology.residues])
+ #print(' Atoms: %s' % [res for res in mytopology.atoms])
+ #print(' Bonds: %s' % [res for res in mytopology.bonds])
+ #print(mytopology.chain(0))
+ table, bonds = mytopology.to_dataframe()
+ print(bonds)
+ print(len(bonds))
+ array = table.get("resSeq")
+ #print([res for res in array])
+ mydict = table.to_dict(orient='list')
+ mylist = mydict["resSeq"]
+ print('mylist')
+ print(mylist)
+ print('mylist')
+ atomindex = np.arange(len(mylist))
+ passarray = np.zeros((len(mylist), 2), dtype=int)
+ passarray[:,0] = atomindex
+ passarray[:,1] = mylist
+ print(table)
+ print(passarray)
+ #print(bonds)
+ #for bond in bonds:
+ # print(bond)
+ mol_to_mol_bond = []
+ mol_to_mol_bond_num = []
+ atom_in_mol_bond = []
+ atom_in_mol_bond_num = []
+ index=0
+ print(len([res for res in mytopology.bonds]))
+ for res in mytopology.bonds:
+ molname1, molatom1 = str(res[0]).split('-')
+ molname2, molatom2 = str(res[1]).split('-')
+ if molname1 in molname2:
+ atom_in_mol_bond.append(res)
+ atom_in_mol_bond_num.append(bonds[index])
+ else:
+ mol_to_mol_bond.append(res)
+ mol_to_mol_bond_num.append(bonds[index])
+ index += 1
+ print(mol_to_mol_bond)
+ print(np.array(mol_to_mol_bond_num))
+
+ atom_list = list(set(mydict["name"]))
+ atom_type_dict = {}
+ print(atom_list)
+ print(len(atom_list))
+ masses = {}
+ for ielm in range(len(atom_list)):
+ elm = atom_list[ielm]
+ atom_type_dict.update({elm : ielm})
+ for atom in mytopology.atoms:
+ if elm == atom.name:
+ masses.update({atom.name : atom.element})
+ print([atom_list, masses])
+ interDict = {}
+ interTypeDict = {}
+ for elm in atom_list:
+ bondList = []
+ typeList = []
+ bondid = 0
+ for bond in mytopology.bonds:
+ molname1, molatom1 = str(bond[0]).split('-')
+ molname2, molatom2 = str(bond[1]).split('-')
+ if (elm == str(molatom1) or elm == str(molatom2)):
+ bondList.append(list(bonds[bondid]))
+ #typeList.append(list([molatom1,molatom2]))
+ typeList.append(list([
+ atom_type_dict[str(molatom1)],
+ atom_type_dict[str(molatom2)]
+ ]))
+ interDict.update({elm : bondList})
+ interTypeDict.update({elm : typeList})
+ bondid += 1
+
+ # print(interDict)
+ print(interDict["CH3"])
+ print(interTypeDict["CH3"])
+ print(len(interDict["CH3"][0]))
+
+ # print([[atom.element.radius, atom.index, atom.name, atom.residue.name, atom.element] for atom in mytopology.atoms])
+
+ # Atom labels:
+ print(table)
+ # for atom in mytopology.atoms:
+ # print(atom)
+ labellist = mydict["element"]
+ #print(labellist)
+ # atomlabels = np.arange(len(labellist))
+ # labelarray = np.zeros((len(labellist), 2), dtype=int)
+ # labelarray[:,0] = atomlabels
+ # labelarray[:,1] = labellist
+ # print(labelarray)
+ # print(readdata.unitcell_lengths)
+ latticevectors = readdata.iread().unitcell_vectors
+ box = readdata.iread().unitcell_lengths
+ step = 0
+ for latvec in latticevectors:
+ print(step, latvec)
+ step += 1
+ for pbc in box:
+ print(step, box)
+ step += 1
+