Ligand charge differs between BioSimSpace and FESetup about biosimspace HOT 7 CLOSED

# formal_charge.py

import openbabel as ob
import sys

# Get the name of the PDB file.
file = sys.argv[1]

# Read the ligand PDB into an OpenBabel molecule.
# The PDB file has no formal charge information.
conv = ob.OBConversion()
mol = ob.OBMol()
conv.ReadFile(mol, file)

# Define some variables from FESetup's "const.py" file.
_phospho = '[#15D4](~[OD1])(~[OD1])(~[OD2])~[OD2]'
_sulfonyl = '[#16](~[OD1])(~[OD1])'

# Again, below is taken from "const.py"
# FIXME: this compensates for missing charge determination by Openbabel, so
#        will fail to detect charges in cases not covered here'; relies on
#        explicit hydrogens!
CHARGE_SMARTS = (
    ('[#16;$([#16D4](~[OD1])(~[OD1])~[OD1])]', -1),  # sulfates, Sac
    ('[#16;$([#16D3](~[OD1])~[OD1])]', -1),  # sulfonates, Sac
    ('[#15;$([#15D4](~[OD1])(~[OD1])~[OD1])]', -2),  # phosphates, Pac
    ('[#15;$([#15D3](~[OD1])~[OD1])]', -2),  # phosphonates, Pac
    ('[#15;$(%s)]' % _phospho, -1),     # polyphosphates or phosphate diesters
    ('[CD3](~[NH2])(~[NH2])~[NH]', +1), # guanidinium
    ('O=C~[CD3,ND2]~C=O', -1),          # beta-dicarbonyls, imides
    ('%s~[ND2]~C(=O)~[NH]' % _sulfonyl, -1), # sulfonylureas
    ('%s~[ND2]~C(=O)~[#6]' % _sulfonyl, -1), # (sulfon)imides
    ('%s~[ND2]~c' % _sulfonyl, -1),     # N-aryl sulfonamides
    ('{0}~[ND2]~{0}'.format(_sulfonyl), -1), # sulfonimides 2
    ('[#6][ND3]([#6])~[#6H]~[ND3]([#6])[#6]', +1), # formamidiniums
    ('C(~[OD1])~[OD1]', -1),            # carbonates, Cac in atomtyp.txt
    ('[#6][#16D1]', -1),                # thiolates
    ('[#6][#8D1]', -1),                 # alkoxides (alcoholates)
    ('[#7D4]', +1),                     # atomtyp.txt checks for [#7X4]
    ('{0}~[ND3](~{0})~[CD3](~{0})~{0}'.format('[#6,#1]'), +1),  # imines
    )

# Zero the total formal charge.
formal_charge = 0

# Initialise the smarts pattern object.
sm = ob.OBSmartsPattern()

# Loop over all smarts strings / formal charge pairs. 
for smarts, fchg in CHARGE_SMARTS:
    sm.Init(smarts)
    
    # Does the molecule match this string.
    # If so, update the formal charge.
    if sm.Match(mol):
        m = list(sm.GetUMapList() )
        formal_charge += fchg * len(m)

# Print the total formal charge.
print(file, formal_charge)