Fix reading AlphaFold files (missing auth_seq_id error)

include/cif++/Structure.hpp

@@ -196,12 +196,18 @@ class Residue
 {
   public:
 	// constructors should be private, but that's not possible for now (needed in emplace)
+
 	// constructor for waters
 	Residue(const Structure& structure, const std::string& compoundID,
 		const std::string& asymID, const std::string& authSeqID);
 
+	// constructor for a residue without a sequence number
+	Residue(const Structure& structure, const std::string& compoundID,
+		const std::string& asymID);
+
+	// constructor for a residue with a sequence number
 	Residue(const Structure& structure, const std::string& compoundID,
-		const std::string& asymID, int seqID = 0);
+		const std::string& asymID, int seqID, const std::string& authSeqID);
 
 	Residue(const Residue& rhs) = delete;
 	Residue& operator=(const Residue& rhs) = delete;
@@ -271,6 +277,9 @@ class Residue
 	const Structure* mStructure = nullptr;
 	std::string	mCompoundID, mAsymID;
 	int mSeqID = 0;
+
+	// Watch out, this is used only to label waters... The rest of the code relies on
+	// MapLabelToAuth to get this info. Perhaps we should rename this member field.
 	std::string mAuthSeqID;	
 	AtomView mAtoms;
 };
@@ -288,8 +297,7 @@ class Monomer : public Residue
 	Monomer(Monomer&& rhs);
 	Monomer& operator=(Monomer&& rhs);
 
-//	Monomer(const Polymer& polymer, uint32_t index);
-	Monomer(const Polymer& polymer, uint32_t index, int seqID,
+	Monomer(const Polymer& polymer, uint32_t index, int seqID, const std::string& authSeqID,
 		const std::string& compoundID);
 
 	bool is_first_in_chain() const;

include/cif++/Symmetry.hpp

@@ -76,6 +76,27 @@ struct SymopData
 		return m_packed == rhs.m_packed;
 	}
 
+	std::array<int,15> data() const
+	{
+		return {
+			static_cast<int>(m_packed >> 34) bitand 0x03,
+			static_cast<int>(m_packed >> 32) bitand 0x03,
+			static_cast<int>(m_packed >> 30) bitand 0x03,
+			static_cast<int>(m_packed >> 28) bitand 0x03,
+			static_cast<int>(m_packed >> 26) bitand 0x03,
+			static_cast<int>(m_packed >> 24) bitand 0x03,
+			static_cast<int>(m_packed >> 22) bitand 0x03,
+			static_cast<int>(m_packed >> 20) bitand 0x03,
+			static_cast<int>(m_packed >> 18) bitand 0x03,
+			static_cast<int>(m_packed >> 15) bitand 0x07,
+			static_cast<int>(m_packed >> 12) bitand 0x07,
+			static_cast<int>(m_packed >>  9) bitand 0x07,
+			static_cast<int>(m_packed >>  6) bitand 0x07,
+			static_cast<int>(m_packed >>  3) bitand 0x07,
+			static_cast<int>(m_packed >>  0) bitand 0x07,
+		};
+	}
+
   private:
 
 	friend struct SymopDataBlock;

src/Cif++.cpp

@@ -2572,8 +2572,6 @@ void Category::update_value(RowSet &&rows, const std::string &tag, const std::st
 						check = std::move(check) && Key(ck) == parent[pk].c_str();
 				}
 
-std::cerr << check << std::endl;
-
 				if (childCat->exists(std::move(check)))	// phew..., narrow escape
 					continue;
 

src/Structure.cpp

@@ -827,12 +827,18 @@ Residue::Residue(const Structure &structure, const std::string &compoundID,
 	assert(not mAtoms.empty());
 }
 
+Residue::Residue(const Structure& structure, const std::string& compoundID, const std::string& asymID)
+	: Residue(structure, compoundID, asymID, 0, {})
+{
+}
+
 Residue::Residue(const Structure &structure, const std::string &compoundID,
-	const std::string &asymID, int seqID)
+	const std::string &asymID, int seqID, const std::string& authSeqID)
 	: mStructure(&structure)
 	, mCompoundID(compoundID)
 	, mAsymID(asymID)
 	, mSeqID(seqID)
+	, mAuthSeqID(authSeqID)
 {
 	assert(mCompoundID != "HOH");
 
@@ -1140,8 +1146,8 @@ std::ostream &operator<<(std::ostream &os, const Residue &res)
 //{
 //}
 
-Monomer::Monomer(const Polymer &polymer, uint32_t index, int seqID, const std::string &compoundID)
-	: Residue(*polymer.structure(), compoundID, polymer.asymID(), seqID)
+Monomer::Monomer(const Polymer &polymer, uint32_t index, int seqID, const std::string& authSeqID, const std::string &compoundID)
+	: Residue(*polymer.structure(), compoundID, polymer.asymID(), seqID, authSeqID)
 	, mPolymer(&polymer)
 	, mIndex(index)
 {
@@ -1640,8 +1646,8 @@ Polymer::Polymer(const Structure &s, const std::string &entityID, const std::str
 	for (auto r : mPolySeq)
 	{
 		int seqID;
-		std::string compoundID;
-		cif::tie(seqID, compoundID) = r.get("seq_id", "mon_id");
+		std::string compoundID, authSeqID;
+		cif::tie(seqID, authSeqID, compoundID) = r.get("seq_id", "auth_seq_num", "mon_id");
 
 		uint32_t index = size();
 
@@ -1649,11 +1655,11 @@ Polymer::Polymer(const Structure &s, const std::string &entityID, const std::str
 		if (not ix.count(seqID))
 		{
 			ix[seqID] = index;
-			emplace_back(*this, index, seqID, compoundID);
+			emplace_back(*this, index, seqID, authSeqID, compoundID);
 		}
 		else if (cif::VERBOSE)
 		{
-			Monomer m{*this, index, seqID, compoundID};
+			Monomer m{*this, index, seqID, authSeqID, compoundID};
 			std::cerr << "Dropping alternate residue " << m << std::endl;
 		}
 	}
@@ -1992,14 +1998,24 @@ std::tuple<char, int, char> Structure::MapLabelToAuth(
 			 cif::Key("asym_id") == asymID and
 			 cif::Key("seq_id") == seqID))
 	{
-		std::string auth_asym_id, pdb_ins_code;
-		int pdb_seq_num;
+		std::string auth_asym_id, pdb_ins_code, pdb_seq_num, auth_seq_num;
+
+		cif::tie(auth_asym_id, auth_seq_num, pdb_seq_num, pdb_ins_code) =
+			r.get("pdb_strand_id", "auth_seq_num", "pdb_seq_num", "pdb_ins_code");
 
-		cif::tie(auth_asym_id, pdb_seq_num, pdb_ins_code) =
-			r.get("pdb_strand_id", "pdb_seq_num", "pdb_ins_code");
+		if (auth_seq_num.empty())
+			auth_seq_num = pdb_seq_num;
 
-		result = std::make_tuple(auth_asym_id.front(), pdb_seq_num,
+		try
+		{
+			result = std::make_tuple(auth_asym_id.front(), std::stoi(auth_seq_num),
+				pdb_ins_code.empty() ? ' ' : pdb_ins_code.front());
+		}
+		catch (const std::exception& ex)
+		{
+			result = std::make_tuple(auth_asym_id.front(), 0,
 				pdb_ins_code.empty() ? ' ' : pdb_ins_code.front());
+		}
 
 		found = true;
 		break;