openbabel/openbabel — pull request #23 by patrickfuller (2013-09-03 03:55:28+00:00)
Title: Improvements to `OBUnitCell::FillUnitCell` and CIF writing

I've been trying to move all my crystallography file i/o to openbabel, and, in its use, have discovered a few bugs. This pull request aims to fix some of them:

**Issue 1:** Density calculations overestimate real value. To replicate, use [CoMOF74](https://gist.github.com/patrickfuller/6419438):

```python
import pybel
with open("CoMOF74.cif") as in_file:
mof = in_file.read()
mol = pybel.readstring("cif", mof)
mol.unitcell.FillUnitCell(mol.OBMol)
# Density in g/cm^3. Should be ~1.2
print sum(atom.atomicmass for atom in mol.atoms) / (mol.unitcell.GetCellVolume() * 0.6022)
```

The problem here is that the duplicate-atom checking in `OBUnitCell::FillUnitCell` assumes the input atoms are already unique. If this isn't the case, the method compounds the error. Fixed by checking the input atoms before expanding the unit cell in `FillUnitCell` (maybe this kind of check belongs somewhere else?).

**Issue 2:** `OBUnitCell::FillUnitCell` hangs on large unit cells.

This is due to the algorithm used in atom checking. It checks atom-atom distances, which is O(N^2). The new approach uses a hash, which should be amortized O(N). It uses the extra speed to check *every* atom for uniqueness, rather than just symmetry unique ones. Note that a better approach would use an octree, but the hash was easier to implement and proves the point.

**Issue 3:** `OBUnitCell::FillUnitCell` doesn't check all possible duplicates. To replicate, use [FeMOF74](
https://gist.github.com/patrickfuller/6419438) and the above script.

The underlying issue is that `FillUnitCell` assumes that the listed symmetry operations are logical. However, CIFs are a mess. In the case of the example cif, the file already has more than one unitcell worth of atoms, and additionally has defined symmetry operations. This is solved by the solution of Issue 2; by using the extra speed afforded by the hash to check every atom-atom connection, we can check for more duplicates.

**Issue 4:** Some crystallographic programs fail if atom coordinates are not in [0, 1].

It's annoying, but an easy fix. I just added `WrapFractionalCoordinate` to cifformat.cpp.

**Issue 5:** Centering can segfault. To replicate, use CoMOF74.

```python
import pybel
with open("CoMOF74.cif") as in_file:
mof = in_file.read()
mol = pybel.readstring("cif", mof)
mol.unitcell.FillUnitCell(mol.OBMol)
# Cycling fixes the issue
# mol = pybel.readstring("cif", mol.write("cif"))
mol.OBMol.Center()
mol.write("cif")
```

The `Center()` method operates off of a bunch of private variables. I'm guessing something isn't getting updated when `FillUnitCell` is called, but I couldn't figure out what it was. I'll look into this one again, but I'd appreciate any pointers in the mean time!