pub_soft/pymol-open-source
3
0
Fork 0
mirror of https://github.com/schrodinger/pymol-open-source.git synced 2026-06-03 19:54:24 +08:00
Code Issues Packages Projects Releases Wiki Activity
Files
master
pymol-open-source/layer2/CoordSet.cpp
Jarrett Johnson 46530027e6 Fix cell color for ray/png 
Fixes #505
2026-05-04 10:58:05 -04:00

1748 lines
49 KiB
C++
Raw Permalink Blame History

/*
A* -------------------------------------------------------------------
B* This file contains source code for the PyMOL computer program
C* copyright 1998-2000 by Warren Lyford Delano of DeLano Scientific.
D* -------------------------------------------------------------------
E* It is unlawful to modify or remove this copyright notice.
F* -------------------------------------------------------------------
G* Please see the accompanying LICENSE file for further information.
H* -------------------------------------------------------------------
I* Additional authors of this source file include:
-*
-*
-*
Z* -------------------------------------------------------------------
*/
#include"os_python.h"
#include"os_numpy.h"
#include"os_std.h"
#include <algorithm>
#include"Base.h"
#include"MemoryDebug.h"
#include"Err.h"
#include"Scene.h"
#include"CoordSet.h"
#include"Color.h"
#include"PConv.h"
#include"P.h"
#include"Matrix.h"
#include"Sphere.h"
#include"Util.h"
#include"Feedback.h"
#include"RepWireBond.h"
#include"RepCylBond.h"
#include"RepDot.h"
#include"RepMesh.h"
#include"RepSphere.h"
#include"RepSphereImmediate.h"
#include"RepRibbon.h"
#include"RepCartoon.h"
#include"RepSurface.h"
#include"RepLabel.h"
#include"RepNonbonded.h"
#include"RepNonbondedSphere.h"
#include"RepEllipsoid.h"
#include"Symmetry.h"
#include"PyMOLGlobals.h"
#include"PyMOLObject.h"
#include "Executive.h"
#include "Lex.h"
#include "Property.h"
/**
* Get atom coordinates, taking symmetry operation into account.
*
* @param v_out Buffer for return value, may or may not be used.
* @param inv If true, apply the inverse of symop
* @return Coordinate pointer, or nullptr if symmetry operation was invalid.
*/
float const* CoordSet::coordPtrSym(
int idx, pymol::SymOp const& symop, float* v_out, bool inv) const
{
auto const* v_in = coordPtr(idx);
// default symop evaluates to false
if (!symop) {
return v_in;
}
auto* sym = getSymmetry();
if (!sym || (symop.index && symop.index >= sym->getNSymMat())) {
return nullptr;
}
copy3f(v_in, v_out);
double const* const statemat = getPremultipliedMatrix();
if (statemat) {
transform44d3f(ObjectStateGetInvMatrix(this), v_out, v_out);
}
transform33f3f(sym->Crystal.realToFrac(), v_out, v_out);
if (inv) {
v_out[0] -= symop.x;
v_out[1] -= symop.y;
v_out[2] -= symop.z;
}
if (symop.index) {
auto mat = sym->getSymMat(symop.index);
if (inv) {
inverse_transform44f3f(mat, v_out, v_out);
} else {
transform44f3f(mat, v_out, v_out);
}
}
if (!inv) {
v_out[0] += symop.x;
v_out[1] += symop.y;
v_out[2] += symop.z;
}
transform33f3f(sym->Crystal.fracToReal(), v_out, v_out);
if (statemat) {
transform44d3f(statemat, v_out, v_out);
}
return v_out;
}
/*========================================================================*/
/**
* Get coordinate index for given atom index
*/
int CoordSet::atmToIdx(int atm) const {
if (Obj->DiscreteFlag) {
if (this == Obj->DiscreteCSet[atm])
return Obj->DiscreteAtmToIdx[atm];
return -1;
}
assert(atm < AtmToIdx.size());
return AtmToIdx[atm];
}
void CoordSet::updateNonDiscreteAtmToIdx(unsigned natom)
{
assert(!Obj || natom == Obj->NAtom);
AtmToIdx.resize(natom);
std::fill_n(AtmToIdx.data(), natom, -1);
for (unsigned idx = 0, idx_end = getNIndex(); idx != idx_end; ++idx) {
auto const atm = IdxToAtm[idx];
assert(atm < natom);
AtmToIdx[atm] = idx;
}
}
/*========================================================================*/
static char sATOM[] = "ATOM ";
static char sHETATM[] = "HETATM";
int CoordSetValidateRefPos(CoordSet * I)
{
if(I->RefPos) {
VLACheck(I->RefPos, RefPosType, I->NIndex);
return true;
} else {
int ok = true && (I->RefPos = pymol::vla<RefPosType>(I->NIndex));
if(ok) {
int a;
for(a = 0; a < I->NIndex; a++) {
const float* src = I->coordPtr(a);
copy3f(src, I->RefPos[a].coord);
I->RefPos[a].specified = true;
}
}
return ok;
}
}
namespace
{
/// Return negative if lhs<rhs, zero if lhs==rhs, positive if lhs>rhs.
inline int cmp(unsigned lhs, unsigned rhs)
{
return lhs < rhs ? -1 : (lhs == rhs ? 0 : 1);
}
inline int cmp(pymol::SymOp const& lhs, pymol::SymOp const& rhs)
{
return cmp( //
reinterpret_cast<unsigned const&>(lhs),
reinterpret_cast<unsigned const&>(rhs));
}
} // namespace
/*========================================================================*/
/**
* Compare bonds by atom indices and symmetry operations. Swapped indices are
* not considered as equal. Bond order, ids and settings are not considered.
*/
int BondCompare(BondType const* a, BondType const* b)
{
int c;
(c = cmp(a->index[0], b->index[0])) == 0 &&
(c = cmp(a->index[1], b->index[1])) == 0 &&
(c = cmp(a->symop_2, b->symop_2));
return c;
}
/*========================================================================*/
int BondInOrder(BondType const* a, int b1, int b2)
{
return (BondCompare(a + b1, a + b2) <= 0);
}
int CoordSetFromPyList(PyMOLGlobals * G, PyObject * list, CoordSet ** cs)
{
CoordSet *I = nullptr;
int ok = true;
int ll = 0;
if(*cs) {
delete *cs;
*cs = nullptr;
}
if(list == Py_None) { /* allow None for CSet */
*cs = nullptr;
} else {
if(ok)
I = CoordSetNew(G);
if(ok)
ok = (I != nullptr);
if(ok)
ok = (list != nullptr);
if(ok)
ok = PyList_Check(list);
if(ok)
ll = PyList_Size(list);
/* TO SUPPORT BACKWARDS COMPATIBILITY...
Always check ll when adding new PyList_GetItem's */
if(ok)
ok = PConvPyIntToInt(PyList_GetItem(list, 0), &I->NIndex);
if(ok)
ok = PConvPyListToFloatVLA(PyList_GetItem(list, 2), &I->Coord);
if(ok){
PConvFromPyListItem(G, list, 3, I->IdxToAtm);
}
if(ok && (ll > 5))
ok = CPythonVal_PConvPyStrToStr_From_List(G, list, 5, I->Name, sizeof(WordType));
if(ok && (ll > 6)){
CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 6);
ok = ObjectStateFromPyList(G, val, I);
CPythonVal_Free(val);
}
if(ok && (ll > 7)){
CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 7);
I->Setting.reset(SettingNewFromPyList(G, val));
CPythonVal_Free(val);
}
if (ok && (ll > 9)) {
CPythonVal* val = CPythonVal_PyList_GetItem(G, list, 9);
I->prop_id = PropertyFromPyList(G, val);
CPythonVal_Free(val);
}
if(ok && (ll > 10)){
CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 10);
if (!CPythonVal_IsNone(val))
I->SculptCGO = CGONewFromPyList(G, val, 0, 1);
else {
I->SculptShaderCGO = I->SculptCGO = nullptr;
}
CPythonVal_Free(val);
}
if(ok){
if(ll > 11){
// load in atom-state settings
CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 11);
if (!CPythonVal_IsNone(val)){
int a;
I->atom_state_setting_id = pymol::vla<int>(I->NIndex);
for (a=0; a<I->NIndex; a++){
CPythonVal *val2 = CPythonVal_PyList_GetItem(G, val, a);
if (!CPythonVal_IsNone(val2)){
CPythonVal_PConvPyIntToInt(val2, &I->atom_state_setting_id[a]);
if (I->atom_state_setting_id[a]) {
I->atom_state_setting_id[a] = SettingUniqueConvertOldSessionID(G, I->atom_state_setting_id[a]);
}
}
CPythonVal_Free(val2);
}
} else {
I->atom_state_setting_id = nullptr;
}
CPythonVal_Free(val);
} else {
// translate legacy LabPos to label_placement_offset
auto val = CPythonVal_PyList_GetItem(G, list, 8);
auto res = CPythonVal_PConvPyListToLabPosVec(G, val);
if (res && !res->empty()) {
auto const& LabPos = *res;
for (int a = 0; a < I->NIndex; ++a) {
if (length3f(LabPos[a].offset) > R_SMALL4) {
SettingSet(
cSetting_label_placement_offset, LabPos[a].offset, I, a);
}
}
}
CPythonVal_Free(val);
}
}
if (ll > 12) {
CPythonVal* val = CPythonVal_PyList_GetItem(G, list, 12);
I->Symmetry.reset(SymmetryNewFromPyList(G, val));
CPythonVal_Free(val);
}
if(!ok) {
delete I;
*cs = nullptr;
} else {
*cs = I;
}
}
return (ok);
}
/**
* Coord set as numpy array
*/
PyObject *CoordSetAsNumPyArray(CoordSet * cs, short copy)
{
#ifndef _PYMOL_NUMPY
PRINTFB(cs->G, FB_CoordSet, FB_Errors)
"No numpy support\n" ENDFB(cs->G);
return nullptr;
#else
PyObject *result = nullptr;
const int base_size = sizeof(float);
int typenum = -1;
npy_intp dims[2] = {0, 3};
import_array1(nullptr);
switch(base_size) {
case 4: typenum = NPY_FLOAT32; break;
case 8: typenum = NPY_FLOAT64; break;
}
if(typenum == -1) {
printf("error: no typenum for float size %d\n", base_size);
return nullptr;
}
dims[0] = cs->NIndex;
if(copy) {
if((result = PyArray_SimpleNew(2, dims, typenum)))
memcpy(PyArray_DATA((PyArrayObject *)result), cs->Coord, cs->NIndex * 3 * base_size);
} else {
result = PyArray_SimpleNewFromData(2, dims, typenum, cs->Coord.data());
}
return result;
#endif
}
PyObject *CoordSetAsPyList(CoordSet * I)
{
PyObject *result = nullptr;
if(I) {
auto G = I->G;
int pse_export_version = SettingGet<float>(G, cSetting_pse_export_version) * 1000;
bool dump_binary = SettingGet<bool>(G, cSetting_pse_binary_dump) && (!pse_export_version || pse_export_version >= 1765);
result = PyList_New(13);
PyList_SetItem(result, 0, PyInt_FromLong(I->NIndex));
int const NAtIndex = I->AtmToIdx.size();
PyList_SetItem(result, 1, PyInt_FromLong(NAtIndex ? NAtIndex : I->Obj->NAtom)); // legacy
PyList_SetItem(result, 2, PConvFloatArrayToPyList(I->Coord, I->NIndex * 3, dump_binary));
PyList_SetItem(result, 3, PConvIntArrayToPyList(I->IdxToAtm.data(), I->NIndex, dump_binary));
if (!I->AtmToIdx.empty()
&& pse_export_version < 1770)
PyList_SetItem(result, 4, PConvIntArrayToPyList(I->AtmToIdx.data(), NAtIndex, dump_binary));
else
PyList_SetItem(result, 4, PConvAutoNone(nullptr));
PyList_SetItem(result, 5, PyString_FromString(I->Name));
PyList_SetItem(result, 6, ObjectStateAsPyList(I));
PyList_SetItem(result, 7, SettingAsPyList(I->Setting.get()));
PyList_SetItem(result, 8, PConvAutoNone(nullptr) /* LabPos */);
PyList_SetItem(result, 9,
I->prop_id ? PropertyAsPyList(G, I->prop_id, true) :
PConvAutoNone(Py_None));
if(I->SculptCGO) {
PyList_SetItem(result, 10, CGOAsPyList(I->SculptCGO));
} else {
PyList_SetItem(result, 10, PConvAutoNone(nullptr));
}
if (I->has_any_atom_state_settings()) {
int a;
PyObject *settings_list = nullptr;
settings_list = PyList_New(I->NIndex);
for (a=0; a<I->NIndex; a++){
if (I->has_atom_state_settings(a)) {
PyList_SetItem(settings_list, a, PyInt_FromLong(I->atom_state_setting_id[a]));
} else {
PyList_SetItem(settings_list, a, PConvAutoNone(nullptr));
}
}
PyList_SetItem(result, 11, settings_list);
} else {
PyList_SetItem(result, 11, PConvAutoNone(nullptr));
}
PyList_SetItem(result, 12, SymmetryAsPyList(I->Symmetry.get()));
/* TODO spheroid, periodic box ... */
}
return (PConvAutoNone(result));
}
/**
* Updates IdxToAtm and adjusts all NIndex sized arrays.
*
* @param lookup atm_old to atm_new index mapping.
* Deleted atoms have `atm_new = -1`.
*
* @pre `lookup[a] <= a`
*/
void CoordSetAdjustAtmIdx(CoordSet* I, const int* lookup)
{
auto G = I->G;
int offset = 0;
for (int idx = 0; idx < I->getNIndex(); ++idx) {
auto const idx_new = idx + offset;
auto const atm_new = lookup[I->IdxToAtm[idx]];
assert(I->IdxToAtm[idx] >= atm_new);
I->IdxToAtm[idx_new] = atm_new;
if (atm_new == -1) {
--offset;
if (I->has_atom_state_settings(idx)) {
SettingUniqueDetachChain(G, I->atom_state_setting_id[idx]);
I->atom_state_setting_id[idx] = 0;
}
} else if (offset) {
copy3f(I->coordPtr(idx), I->coordPtr(idx_new));
if (I->RefPos) {
I->RefPos[idx_new] = std::move(I->RefPos[idx]);
}
if (I->has_atom_state_settings(idx)) {
I->atom_state_setting_id[idx_new] = std::move(I->atom_state_setting_id[idx]);
I->atom_state_setting_id[idx] = 0;
}
}
}
assert(offset <= 0);
if (offset < 0) {
I->setNIndex(I->getNIndex() + offset);
I->invalidateRep(cRepAll, cRepInvAtoms); /* this will free Color */
}
}
CoordSet* CoordSetCopy(const CoordSet* src)
{
if(!src) {
return nullptr;
}
return new CoordSet(*src);
}
/*========================================================================*/
/**
* Append cs to I.
*
* @pre I and cs are non-overlapping
*/
int CoordSetMerge(ObjectMolecule *OM, CoordSet * I, const CoordSet * cs)
{
assert(OM == I->Obj);
auto const nIndexOld = I->getNIndex();
I->setNIndex(I->getNIndex() + cs->getNIndex());
for (int idx_src = 0; idx_src < cs->getNIndex(); ++idx_src) {
int const idx = idx_src + nIndexOld;
int const atm = cs->IdxToAtm[idx_src];
I->IdxToAtm[idx] = cs->IdxToAtm[idx_src];
if (OM->DiscreteFlag) {
OM->DiscreteAtmToIdx[atm] = idx;
OM->DiscreteCSet[atm] = I;
} else {
I->AtmToIdx[atm] = idx;
}
copy3f(cs->coordPtr(idx_src), I->coordPtr(idx));
}
if (cs->RefPos) {
I->RefPos.resize(I->getNIndex());
std::copy_n(cs->RefPos.data(), cs->getNIndex(), I->RefPos.data() + nIndexOld);
}
// TODO copy or move atom_state_setting_id here?
I->invalidateRep(cRepAll, cRepInvAll);
return true;
}
/*========================================================================*/
int CoordSetTransformAtomTTTf(CoordSet * I, int at, const float *TTT)
{
int a1 = I->atmToIdx(at);
if(a1 < 0)
return false;
auto* v1 = I->coordPtr(a1);
MatrixTransformTTTfN3f(1, v1, TTT, v1);
return true;
}
/*========================================================================*/
int CoordSetTransformAtomR44f(CoordSet * I, int at, const float *matrix)
{
int a1 = I->atmToIdx(at);
if(a1 < 0)
return false;
auto* v1 = I->coordPtr(a1);
MatrixTransformR44fN3f(1, v1, matrix, v1);
return true;
}
/*========================================================================*/
void CoordSetRecordTxfApplied(CoordSet * I, const float *matrix, int homogenous)
{
double temp[16];
if(!homogenous) {
convertTTTfR44d(matrix, temp);
} else {
convert44f44d(matrix, temp);
}
ObjectStateLeftCombineMatrixR44d(I, temp);
}
/*========================================================================*/
int CoordSetMoveAtom(CoordSet * I, int at, const float *v, int mode)
{
int a1 = I->atmToIdx(at);
if(a1 < 0)
return false;
auto* v1 = I->coordPtr(a1);
if(mode) {
add3f(v, v1, v1);
} else {
copy3f(v, v1);
}
return true;
}
/*========================================================================*/
/**
* Retrieve the label offset for a given atom
* @param I target Coordset
* @param atm atom index
*/
pymol::Result<pymol::Vec3> CoordSet::getAtomLabelOffset(int atm) const
{
int idx = this->atmToIdx(atm);
if (idx < 0) {
return pymol::make_error("Invalid atom Idx");
}
auto G = this->G;
auto obj = this->Obj;
pymol::Vec3 result{};
const float* at_offset_ptr;
auto ai = &obj->AtomInfo[atm];
int at_label_relative_mode = 0;
AtomStateGetSetting_i(G, obj, this, idx, ai, cSetting_label_relative_mode,
&at_label_relative_mode);
switch (at_label_relative_mode) {
case cLabelRelativeMode::Default:
AtomStateGetSetting(
G, obj, this, idx, ai, cSetting_label_placement_offset, &at_offset_ptr);
break;
case cLabelRelativeMode::ScreenRelative:
case cLabelRelativeMode::ScreenPixelSpace:
AtomStateGetSetting(
G, obj, this, idx, ai, cSetting_label_screen_point, &at_offset_ptr);
break;
}
std::copy_n(at_offset_ptr, result.size(), result.data());
return result;
}
/**
* Retrieve the label offset for a given atom
* @param I target Coordset
* @param atm atom index
* @param offset offset to be set at
*/
pymol::Result<> CoordSet::setAtomLabelOffset(
int atm, const float* offset)
{
int idx = this->atmToIdx(atm);
if (idx < 0) {
return pymol::make_error("Invalid atom Idx");
}
auto G = this->G;
auto obj = this->Obj;
auto ai = &obj->AtomInfo[atm];
int at_label_relative_mode = 0;
AtomStateGetSetting_i(G, obj, this, idx, ai, cSetting_label_relative_mode,
&at_label_relative_mode);
switch (at_label_relative_mode) {
case cLabelRelativeMode::Default:
SettingSet(cSetting_label_placement_offset, offset, this, idx);
case cLabelRelativeMode::ScreenRelative:
case cLabelRelativeMode::ScreenPixelSpace:
SettingSet(cSetting_label_screen_point, offset, this, idx);
break;
}
return {};
}
void CoordSet::setTitle(pymol::zstring_view title)
{
UtilNCopy(Name, title.c_str(), sizeof(WordType));
}
/**
* Get the transformation matrix which is pre-multiplied to the coordiantes, or
* nullptr if there is either no matrix set, or it's not pre-multiplied
* (matrix_mode=1).
*/
double const* CoordSet::getPremultipliedMatrix() const {
return (SettingGet<int>(*this, cSetting_matrix_mode) > 0)
? nullptr
: ObjectStateGetMatrix(this);
}
int CoordSetMoveAtomLabel(CoordSet * I, int at, const float *v, const float *diff)
{
auto G = I->G;
ObjectMolecule *obj = I->Obj;
int a1 = I->atmToIdx(at);
int result = 0;
/* if label is valid, get the label offset
* and set the new position relative to that */
if(a1 >= 0) {
float at_offset[3];
const float * at_offset_ptr;
int at_label_relative_mode = 0;
AtomInfoType *ai = obj->AtomInfo + at;
AtomStateGetSetting_i(G, obj, I, a1, ai, cSetting_label_relative_mode, &at_label_relative_mode);
switch (at_label_relative_mode){
case cLabelRelativeMode::Default:
AtomStateGetSetting(G, obj, I, a1, ai, cSetting_label_placement_offset, &at_offset_ptr);
add3f(v, at_offset_ptr, at_offset);
SettingSet(cSetting_label_placement_offset, at_offset, I, a1);
break;
case cLabelRelativeMode::ScreenRelative: // screen relative
case cLabelRelativeMode::ScreenPixelSpace: // screen pixel space
{
float voff[3];
int width, height;
SceneGetWidthHeight(G, &width, &height);
if (at_label_relative_mode==1){
voff[0] = 2.f * diff[0] / width;
voff[1] = 2.f * diff[1] / height;
} else {
voff[0] = diff[0];
voff[1] = diff[1];
}
voff[2] = 0.f;
AtomStateGetSetting(G, obj, I, a1, ai, cSetting_label_screen_point, &at_offset_ptr);
add3f(voff, at_offset_ptr, at_offset);
SettingSet(cSetting_label_screen_point, at_offset, I, a1);
}
break;
}
}
return (result);
}
/*========================================================================*/
int CoordSetGetAtomVertex(const CoordSet * I, int at, float *v)
{
int a1 = I->atmToIdx(at);
if(a1 < 0)
return false;
copy3f(I->coordPtr(a1), v);
return true;
}
/*========================================================================*/
int CoordSetGetAtomTxfVertex(const CoordSet * I, int at, float *v)
{
ObjectMolecule *obj = I->Obj;
int a1 = I->atmToIdx(at);
if(a1 < 0)
return false;
copy3f(I->coordPtr(a1), v);
/* apply state transformation */
if (!I->Matrix.empty() && SettingGet<int>(*I, cSetting_matrix_mode) > 0) {
transform44d3f(I->Matrix.data(), v, v);
}
/* object transformation */
if(obj->TTTFlag) {
transformTTT44f3f(obj->TTT, v, v);
}
return true;
}
/*========================================================================*/
int CoordSetSetAtomVertex(CoordSet * I, int at, const float *v)
{
int a1 = I->atmToIdx(at);
if(a1 < 0)
return false;
copy3f(v, I->coordPtr(a1));
return true;
}
/*========================================================================*/
void CoordSetRealToFrac(CoordSet * I, const CCrystal * cryst)
{
if (I->getPremultipliedMatrix()) {
float mat[16];
copy44d44f(ObjectStateGetInvMatrix(I), mat);
CoordSetTransform44f(I, mat);
}
CoordSetTransform33f(I, cryst->realToFrac());
}
/*========================================================================*/
void CoordSetTransform44f(CoordSet * I, const float *mat)
{
int a;
float* v = I->Coord.data();
for(a = 0; a < I->NIndex; a++) {
transform44f3f(mat, v, v);
v += 3;
}
}
/*========================================================================*/
void CoordSetTransform33f(CoordSet * I, const float *mat)
{
int a;
float* v = I->Coord.data();
for(a = 0; a < I->NIndex; a++) {
transform33f3f(mat, v, v);
v += 3;
}
}
/*========================================================================*/
void CoordSetGetAverage(const CoordSet * I, float *v0)
{
int a;
double accum[3];
if(I->NIndex) {
const float* v = I->Coord.data();
accum[0] = *(v++);
accum[1] = *(v++);
accum[2] = *(v++);
for(a = 1; a < I->NIndex; a++) {
accum[0] += *(v++);
accum[1] += *(v++);
accum[2] += *(v++);
}
v0[0] = (float) (accum[0] / I->NIndex);
v0[1] = (float) (accum[1] / I->NIndex);
v0[2] = (float) (accum[2] / I->NIndex);
}
}
/*========================================================================*/
void CoordSetFracToReal(CoordSet * I, const CCrystal * cryst)
{
CoordSetTransform33f(I, cryst->fracToReal());
if (double const* statemat = I->getPremultipliedMatrix()) {
float mat[16];
copy44d44f(statemat, mat);
CoordSetTransform44f(I, mat);
}
}
/**
* Apply the `sca` (SCALEn) transformation to transform to fractional space,
* and then the crystals `FracToReal` transformation to transform back to
* cartesian space.
*
* Don't do anything if pdb_insure_orthogonal=off.
*
* Don't do anything if SCALEn or CRYST1 look bogus. There is a number of
* structures in the PDB which have meaningless values for those.
*
* Without this, creating symmetry mates might produce wrong results.
*/
bool CoordSetInsureOrthogonal(PyMOLGlobals * G,
CoordSet * cset, // coord set to modify
const float * sca, // 4x4 SCALE
const CCrystal *cryst,
bool quiet)
{
if (!SettingGet<bool>(G, cSetting_pdb_insure_orthogonal))
return false;
if (!cryst)
cryst = &cset->Symmetry->Crystal;
auto const r2f = cryst->realToFrac();
// are the matrices sufficiently close to be the same?
if (!sca[3] && !sca[7] && !sca[11] &&
is_allclosef(3, r2f, 3, sca, 4, R_SMALL4)) {
return false;
}
// is the cell a orthogonal 1x1x1? If so, then it should probably be ignored...
// is SCALEn the identity matrix? If so, then it should probably be ignored...
if (is_identityf(3, r2f, R_SMALL4) ||
is_identityf(4, sca, R_SMALL4)) {
PRINTFB(G, FB_ObjectMolecule, FB_Blather)
" ObjectMolReadPDBStr: ignoring SCALEn (identity matrix).\n" ENDFB(G);
return false;
}
// is SCALEn invalid? If so, then it should definitely be ignored...
if (determinant33f(sca, 4) < R_SMALL8 ||
determinant33f(r2f, 3) < R_SMALL8) {
PRINTFB(G, FB_ObjectMolecule, FB_Blather)
" ObjectMolReadPDBStr: ignoring SCALEn (invalid matrix).\n" ENDFB(G);
return false;
}
PRINTFB(G, FB_ObjectMolecule, quiet ? FB_Blather : FB_Actions)
" ObjectMolecule: using SCALEn to compute orthogonal coordinates.\n"
ENDFB(G);
CoordSetTransform44f(cset, sca);
CoordSetFracToReal(cset, cryst);
return true;
}
/*========================================================================*/
/* Rotates the ANISOU vector
*
* matrix: flat 4x4, but only rotation (upper left 3x3) is considered
* anisou: has 6 elements (of symmetric 3x3) and will be rotated in-place
*/
bool RotateU(const double *matrix, float *anisou)
{
int i, j, k;
float Re[3][3];
double e_val[3], e_vec[3][3];
double U[3][3] = {
{ anisou[0], anisou[3], anisou[4] },
{ anisou[3], anisou[1], anisou[5] },
{ anisou[4], anisou[5], anisou[2] },
};
// e_val, e_vec = linalg.eigh(U)
if(!xx_matrix_jacobi_solve(*e_vec, e_val, &i, *U, 3))
return false;
// Re = dot(matrix[:3,:3], e_vec)
for (i = 0; i < 3; i++)
for (j = 0; j < 3; j++) {
Re[i][j] = 0.0;
for (k = 0; k < 3; k++)
Re[i][j] += matrix[i * 4 + k] * e_vec[k][j];
}
// U = dot(Re * e_val, Re.T)
for (i = 0; i < 3; i++)
for (j = 0; j <= i; j++) {
U[i][j] = 0.0;
for (k = 0; k < 3; k++)
U[i][j] += Re[i][k] * e_val[k] * Re[j][k];
}
anisou[0] = U[0][0];
anisou[1] = U[1][1];
anisou[2] = U[2][2];
anisou[3] = U[1][0];
anisou[4] = U[2][0];
anisou[5] = U[2][1];
return true;
}
/*========================================================================*/
/**
* @param v 3x1 vertex in final output space
* @param matrix 4x4 homogenous transformation matrix from model space to output
* space (view matrix * state matrix). Used for ANISOU.
*/
void CoordSetAtomToPDBStrVLA(PyMOLGlobals * G, char **charVLA, int *c,
const AtomInfoType * ai,
const float *v, int cnt,
const PDBInfoRec * pdb_info,
const double *matrix)
{
char *aType;
AtomName name;
ResName resn;
lexidx_t chain;
char formalCharge[4];
bool const ignore_pdb_segi = SettingGet<bool>(G, cSetting_ignore_pdb_segi);
WordType x, y, z;
AtomInfoGetAlignedPDBResidueName(G, ai, resn);
AtomInfoGetAlignedPDBAtomName(G, ai, resn, name);
formalCharge[0] = 0;
if (SettingGet<bool>(G, cSetting_pdb_formal_charges)) {
if((ai->formalCharge > 0) && (ai->formalCharge < 10)) {
sprintf(formalCharge, "%d+", ai->formalCharge);
} else if((ai->formalCharge < 0) && (ai->formalCharge > -10)) {
sprintf(formalCharge, "%d-", -ai->formalCharge);
}
}
if(ai->hetatm)
aType = sHETATM;
else
aType = sATOM;
char inscode = ai->getInscode(true);
VLACheck(*charVLA, char, (*c) + 1000);
if (SettingGet<bool>(G, cSetting_pdb_retain_ids)) {
cnt = ai->id - 1;
}
if(cnt > 99998)
cnt = 99998;
if((!pdb_info) || (!pdb_info->is_pqr_file())) { /* relying upon short-circuit */
short linelen;
sprintf(x, "%8.3f", v[0]);
x[8] = 0;
sprintf(y, "%8.3f", v[1]);
y[8] = 0;
sprintf(z, "%8.3f", v[2]);
z[8] = 0;
linelen =
sprintf((*charVLA) + (*c),
"%6s%5i %-4s%1s%-4s%1.1s%4i%c %s%s%s%6.2f%6.2f %-4.4s%2s%2s\n", aType,
cnt + 1, name, ai->alt, resn, LexStr(G, ai->chain), ai->resv % 10000, inscode, x, y, z, ai->q, ai->b,
ignore_pdb_segi ? "" :
LexStr(G, ai->segi), ai->elem, formalCharge);
if(ai->anisou) {
// Columns 7 - 27 and 73 - 80 are identical to the corresponding ATOM/HETATM record.
char *atomline = (*charVLA) + (*c);
char *anisoline = atomline + linelen;
float anisou[6];
std::copy_n(ai->anisou, 6, anisou);
if(matrix && !RotateU(matrix, anisou)) {
PRINTFB(G, FB_CoordSet, FB_Errors) "RotateU failed\n" ENDFB(G);
return;
}
strncpy(anisoline + 6, atomline + 6, 22);
sprintf(anisoline + 28,
"%7.0f%7.0f%7.0f%7.0f%7.0f%7.0f",
anisou[0] * 1e4, anisou[1] * 1e4, anisou[2] * 1e4,
anisou[3] * 1e4, anisou[4] * 1e4, anisou[5] * 1e4);
strcpy(anisoline + 70, atomline + 70);
strncpy(anisoline, "ANISOU", 6);
(*c) += linelen;
}
(*c) += linelen;
} else {
Chain alt;
if(pdb_info->is_pqr_file() && pdb_info->pqr_workarounds) {
inscode = ' ';
chain = 0; /* no chain IDs */
alt[0] = 0; /* not alt conf identifiers */
} else {
alt[0] = ai->alt[0];
alt[1] = 0;
chain = ai->chain;
}
if (pymol::zstring_view(resn).find_first_not_of(' ') ==
pymol::zstring_view::npos) {
// APBS fails with empty residue names
assert(resn[0] == ' ');
resn[0] = '.';
}
sprintf(x, "%8.3f", v[0]);
if(x[0] != 32)
sprintf(x, " %7.2f", v[0]);
x[8] = 0;
sprintf(y, "%8.3f", v[1]);
y[8] = 0;
if(y[0] != 32)
sprintf(y, " %7.2f", v[1]);
y[8] = 0;
sprintf(z, "%8.3f", v[2]);
if(z[0] != 32)
sprintf(z, " %7.2f", v[2]);
z[8] = 0;
(*c) += sprintf((*charVLA) + (*c), "%6s%5i %-4s%1s%-4s%1.1s%4i%c %s%s%s %11.8f %7.3f\n",
aType, cnt + 1, name, alt, resn,
LexStr(G, chain), ai->resv, inscode, x, y, z, ai->partialCharge, ai->elec_radius);
}
}
/*========================================================================*/
PyObject *CoordSetAtomToChemPyAtom(PyMOLGlobals * G, AtomInfoType * ai, ObjectMolecule *obj,
const float *v, const float *ref, int index, const double *matrix)
{
#ifdef _PYMOL_NOPY
return nullptr;
#else
PyObject *atom = PYOBJECT_CALLMETHOD(P_chempy, "Atom", "");
if(!atom)
ErrMessage(G, "CoordSetAtomToChemPyAtom", "can't create atom");
else {
float tmp_array[6] = { 0.f, 0.f, 0.f, 0.f, 0.f, 0.f };
if (ai->anisou) {
std::copy_n(ai->anisou, 6, tmp_array);
if (matrix)
RotateU(matrix, tmp_array);
}
PConvStringToPyObjAttr(atom, "model", obj->Name);
PConvFloat3ToPyObjAttr(atom, "coord", v);
if(ref)
PConvFloat3ToPyObjAttr(atom, "ref_coord", ref);
if (ai->name)
PConvStringToPyObjAttr(atom, "name", LexStr(G, ai->name));
PConvStringToPyObjAttr(atom, "symbol", ai->elem);
// TODO defaults to UNK if (ai->resn)
PConvStringToPyObjAttr(atom, "resn", LexStr(G, ai->resn));
if (ai->inscode) {
char ins_code[2] = {ai->inscode, 0};
PConvStringToPyObjAttr(atom, "ins_code", ins_code);
}
if (ai->ssType[0])
PConvStringToPyObjAttr(atom, "ss", ai->ssType);
// TODO defaults to 1 if (ai->resv)
PConvIntToPyObjAttr(atom, "resi_number", ai->resv);
if (ai->stereo)
PConvIntToPyObjAttr(atom, "stereo", ai->stereo);
if (ai->chain)
PConvStringToPyObjAttr(atom, "chain", LexStr(G, ai->chain));
if(ai->alt[0])
PConvStringToPyObjAttr(atom, "alt", ai->alt);
if (ai->segi)
PConvStringToPyObjAttr(atom, "segi", LexStr(G, ai->segi));
if (ai->q != 1.)
PConvFloatToPyObjAttr(atom, "q", ai->q);
if (ai->b)
PConvFloatToPyObjAttr(atom, "b", ai->b);
if (ai->anisou)
{
{
PyObject *tmp_obj = PConvFloatArrayToPyList(tmp_array, 6);
if(tmp_obj) {
PyObject_SetAttrString(atom, "u_aniso", tmp_obj);
Py_XDECREF(tmp_obj);
}
}
}
PConvFloatToPyObjAttr(atom, "vdw", ai->vdw);
if (ai->elec_radius)
PConvFloatToPyObjAttr(atom, "elec_radius", ai->elec_radius);
if (ai->partialCharge)
PConvFloatToPyObjAttr(atom, "partial_charge", ai->partialCharge);
if (ai->formalCharge)
PConvIntToPyObjAttr(atom, "formal_charge", ai->formalCharge);
// TODO customType=0 from most files
if(ai->customType != -9999)
PConvIntToPyObjAttr(atom, "numeric_type", ai->customType);
if (ai->textType)
PConvStringToPyObjAttr(atom, "text_type", LexStr(G, ai->textType));
if (ai->custom)
PConvStringToPyObjAttr(atom, "custom", LexStr(G, ai->custom));
PConvIntToPyObjAttr(atom, "hetatm", ai->hetatm);
PConvIntToPyObjAttr(atom, "flags", ai->flags);
PConvIntToPyObjAttr(atom, "id", ai->id); /* not necc. unique */
PConvIntToPyObjAttr(atom, "index", index + 1); /* fragile */
if (ai->prop_id) {
PyObject* props = PropertyAsPyList(G, ai->prop_id, false);
if (PyList_Check(props)) {
PyObject* atomProperties =
PyObject_GetAttrString(atom, "atom_properties");
int lp, lsz = PyList_Size(props);
for (lp = 0; lp < lsz; lp++) {
PyObject* atProp = PyList_GetItem(props, lp);
PyDict_SetItem(atomProperties, PyList_GetItem(atProp, 0),
PyList_GetItem(atProp, 1));
}
}
}
}
if(PyErr_Occurred())
PyErr_Print();
return (atom);
#endif
}
/*========================================================================*/
void CoordSet::invalidateRep(cRep_t type, cRepInv_t level)
{
if(level >= cRepInvVisib) {
if (Obj)
Obj->RepVisCacheValid = false;
}
/* graphical representations need redrawing */
if(level == cRepInvVisib) {
/* cartoon_side_chain_helper */
if (SettingGet<bool>(*this, cSetting_cartoon_side_chain_helper)) {
if((type == cRepCyl) || (type == cRepLine) || (type == cRepSphere))
invalidateRep(cRepCartoon, cRepInvVisib2);
else if(type == cRepCartoon) {
invalidateRep(cRepLine, cRepInvVisib2);
invalidateRep(cRepCyl, cRepInvVisib2);
invalidateRep(cRepSphere, cRepInvVisib2);
}
}
/* ribbon_side_chain_helper */
if (SettingGet<bool>(*this, cSetting_ribbon_side_chain_helper)) {
if((type == cRepCyl) || (type == cRepLine) || (type == cRepSphere))
invalidateRep(cRepRibbon, cRepInvVisib2);
else if(type == cRepRibbon) {
invalidateRep(cRepLine, cRepInvVisib2);
invalidateRep(cRepCyl, cRepInvVisib2);
invalidateRep(cRepSphere, cRepInvVisib2);
}
}
/* line_stick helper */
if (SettingGet<bool>(*this, cSetting_line_stick_helper)) {
if(type == cRepCyl)
invalidateRep(cRepLine, cRepInvVisib2);
else if(type == cRepLine) {
invalidateRep(cRepCyl, cRepInvVisib2);
}
}
}
// cell
if (type == cRepCell) {
UnitCellCGO.reset();
}
/* invalidate basd on one representation, 'type' */
for (RepIterator iter(G, type); iter.next(); ){
auto eff_level = level;
auto const a = iter.rep;
if(level == cRepInvPick) {
switch (a) {
case cRepSurface:
case cRepMesh:
case cRepDot:
/* skip the expensive to recompute, non-pickable
representations */
break;
default: /* default behavior is to blow away the representation */
eff_level = cRepInvRep;
break;
}
}
if(eff_level >= cRepInvVisib) /* make active if visibility has changed */
Active[a] = true;
if (Rep[a]) {
if (eff_level < cRepInvPurge) {
Rep[a]->invalidate(eff_level);
} else {
delete Rep[a];
Rep[a] = nullptr;
}
}
}
if(level >= cRepInvCoord) { /* if coordinates change, then this map becomes invalid */
MapFree(Coord2Idx);
Coord2Idx = nullptr;
ExecutiveInvalidateSelectionIndicatorsCGO(G);
SceneInvalidatePicking(G);
/* invalidate distances */
}
#ifndef NO_MMLIBS
if (level >= cRepInvProp) {
if (validMMStereo == MMPYMOLX_PROP_STATE_AUTO)
validMMStereo = MMPYMOLX_PROP_STATE_NULL;
if (validTextType == MMPYMOLX_PROP_STATE_AUTO)
validTextType = MMPYMOLX_PROP_STATE_NULL;
}
#endif
SceneChanged(G);
}
/*========================================================================*/
#define RepUpdateMacro(rep, new_fn, state) \
{ \
if (Active[rep] && !G->Interrupt) { \
if (Rep[rep]) { \
assert(Rep[rep]->cs == this); \
assert(Rep[rep]->getState() == state); \
Rep[rep] = Rep[rep]->update(); \
} else { \
Rep[rep] = new_fn(this, state); \
if (Rep[rep]) { \
Rep[rep]->fNew = new_fn; \
SceneInvalidatePicking(G); \
} else { \
Active[rep] = false; \
} \
} \
} \
OrthoBusyFast(G, rep, cRepCnt); \
}
/**
* @brief resolve the unit cell color, falling back to the object color when
* `cell_color` is unset (negative).
* @param csSetting pointer to coordset setting
* @param objSetting pointer to object setting
* @param objColor object Color id
* @return color id used to color cell
*/
static int ResolveCellColor(PyMOLGlobals* G, CSetting const* csSetting,
CSetting const* objSetting, int objColor)
{
auto cell_color =
SettingGet_color(G, csSetting, objSetting, cSetting_cell_color);
if (cell_color < 0) {
cell_color = objColor;
}
return cell_color;
}
/*========================================================================*/
void CoordSet::update(int state)
{
assert(G == Obj->G);
OrthoBusyFast(G, 0, cRepCnt);
RepUpdateMacro(cRepLine, RepWireBondNew, state);
RepUpdateMacro(cRepCyl, RepCylBondNew, state);
RepUpdateMacro(cRepDot, RepDotNew, state);
RepUpdateMacro(cRepMesh, RepMeshNew, state);
RepUpdateMacro(cRepSphere, RepSphereNew, state);
RepUpdateMacro(cRepRibbon, RepRibbonNew, state);
RepUpdateMacro(cRepCartoon, RepCartoonNew, state);
RepUpdateMacro(cRepSurface, RepSurfaceNew, state);
RepUpdateMacro(cRepLabel, RepLabelNew, state);
RepUpdateMacro(cRepNonbonded, RepNonbondedNew, state);
RepUpdateMacro(cRepNonbondedSphere, RepNonbondedSphereNew, state);
RepUpdateMacro(cRepEllipsoid, RepEllipsoidNew, state);
for (int a = 0; a < cRepCnt; ++a) {
if (!Rep[a])
Active[a] = false;
}
// cell
if ((Obj->visRep & cRepCellBit) && !UnitCellCGO) {
if (auto const* sym = getSymmetry()) {
UnitCellCGO.reset(CrystalGetUnitCellCGO(&sym->Crystal));
auto use_shader = SettingGet<bool>(G, cSetting_use_shaders);
if (use_shader) {
auto const cell_color = ResolveCellColor(
G, this->Setting.get(), Obj->Setting.get(), Obj->Color);
auto const* color = ColorGet(G, cell_color);
auto preCGO = std::make_unique<CGO>(G);
CGOColorv(preCGO.get(), color);
CGOAppendNoStop(preCGO.get(), UnitCellCGO.get());
std::unique_ptr<CGO> optimized(CGOOptimizeToVBONotIndexed(preCGO.get(), 0));
UnitCellShaderCGO.reset(optimized.release());
assert(UnitCellShaderCGO->use_shader);
}
}
}
SceneInvalidate(G);
OrthoBusyFast(G, 1, 1);
}
/*========================================================================*/
void CoordSetUpdateCoord2IdxMap(CoordSet * I, float cutoff)
{
if(cutoff < R_SMALL4)
cutoff = R_SMALL4;
if(I->NIndex > 10) {
if(I->Coord2Idx) {
if((I->Coord2IdxDiv < cutoff) ||
(((cutoff - I->Coord2IdxReq) / I->Coord2IdxReq) < -0.5F)) {
MapFree(I->Coord2Idx);
I->Coord2Idx = nullptr;
}
}
if(I->NIndex && (!I->Coord2Idx)) { /* NOTE: map based on stored coords */
I->Coord2IdxReq = cutoff;
I->Coord2IdxDiv = cutoff * 1.25F;
I->Coord2Idx = new MapType(I->G, I->Coord2IdxDiv, I->Coord, I->NIndex, nullptr);
if(I->Coord2IdxDiv < I->Coord2Idx->Div)
I->Coord2IdxDiv = I->Coord2Idx->Div;
}
}
}
/*========================================================================*/
void CoordSet::render(RenderInfo * info)
{
const auto pass = info->pass;
const auto pick = bool(info->pick);
auto* ray = info->ray;
auto use_shader = SettingGet<bool>(G, cSetting_use_shaders);
if (!(ray || pick) &&
(SettingGet<int>(*this, cSetting_defer_builds_mode) == 5)) {
if (pass == RenderPass::Antialias) {
ObjectUseColor(Obj);
if (Active[cRepLine])
RepWireBondRenderImmediate(this, info);
if (Active[cRepNonbonded])
RepNonbondedRenderImmediate(this, info);
if (Active[cRepSphere])
RepSphereRenderImmediate(this, info);
if (Active[cRepCyl])
RepCylBondRenderImmediate(this, info);
if (Active[cRepRibbon])
RepRibbonRenderImmediate(this, info);
}
return;
}
const auto sculpt_vdw_vis_mode =
SettingGet<int>(*this, cSetting_sculpt_vdw_vis_mode);
if (pass == RenderPass::Antialias && sculpt_vdw_vis_mode && SculptCGO &&
(Obj->visRep & cRepCGOBit)) {
if (ray) {
CGORenderRay(SculptCGO, ray, info, ColorGet(G, Obj->Color), nullptr,
Setting.get(), Obj->Setting.get());
} else if (G->HaveGUI && G->ValidContext && !pick) {
if (!info->use_shaders) {
CGOFree(SculptShaderCGO);
} else if (!SculptShaderCGO) {
SculptShaderCGO = CGOOptimizeToVBONotIndexed(SculptCGO);
}
auto* cgo = SculptShaderCGO ? SculptShaderCGO : SculptCGO;
CGORender(
cgo, nullptr, Setting.get(), Obj->Setting.get(), info, nullptr);
}
}
// cell
if (UnitCellCGO && (Obj->visRep & cRepCellBit)) {
auto const cell_color = ResolveCellColor(
G, this->Setting.get(), Obj->Setting.get(), Obj->Color);
auto const* color = ColorGet(G, cell_color);
if (ray) {
CGORenderRay(UnitCellCGO.get(), ray, info, color,
nullptr, Setting.get(), Obj->Setting.get());
} else if (!pick && pass == RenderPass::Opaque && G->HaveGUI &&
G->ValidContext) {
auto renderCGO = UnitCellCGO.get();
if (use_shader && UnitCellShaderCGO) {
renderCGO = UnitCellShaderCGO.get();
}
CGORender(renderCGO, color, Setting.get(),
Obj->Setting.get(), info, nullptr);
}
}
auto aastart = cRep_t(0), aaend = cRepCnt;
if (pick) {
int pick_labels = SettingGet<int>(*this, cSetting_pick_labels);
if (pick_labels == 2) { // only pick labels
aastart = cRepLabel;
aaend = cRep_t(aastart + 1);
}
}
for (auto a = aastart; a != aaend; ++a) {
auto* const r = Rep[a];
if (!(r && Active[a])) {
continue;
}
if (!ray) {
ObjectUseColor(Obj);
} else {
ray->wobble(SettingGet<int>(*this, cSetting_ray_texture),
SettingGet<const float*>(*this, cSetting_ray_texture_settings));
ray->color3fv(ColorGet(G, Obj->Color));
}
if (ray || pick) {
r->render(info);
continue;
}
// OIT transparency mode
bool const t_mode_3 = (3 == SettingGet<int>(G, cSetting_transparency_mode));
// render both opaque and transparent pass
bool render_both = t_mode_3;
switch (a) {
case cRepVolume: {
if (t_mode_3 && pass == RenderPass::Transparent) {
r->render(info);
}
} break;
case cRepLabel:
if (pass == RenderPass::Transparent ||
(t_mode_3 && pass == RenderPass::Opaque /* TODO_OPENVR 0 */)) {
r->render(info);
}
break;
case cRepDot:
case cRepCGO:
case cRepCallback:
if (pass == RenderPass::Opaque) {
r->render(info);
}
break;
case cRepCell:
// not implemented
assert(false);
case cRepLine:
case cRepMesh:
case cRepDash:
case cRepExtent:
if (pass == RenderPass::Antialias) {
r->render(info);
}
break;
case cRepNonbonded:
case cRepRibbon:
render_both = false; // cartoon and nonbonded do not have atom-level
// transparency
case cRepNonbondedSphere:
case cRepSurface:
case cRepEllipsoid:
case cRepCyl:
case cRepCartoon:
case cRepSphere:
if (render_both) {
if (pass != RenderPass::Antialias) {
r->render(info);
}
} else if (info->alpha_cgo && a == cRepSurface) {
if (pass == RenderPass::Opaque) {
r->render(info);
}
} else if (r->hasTransparency()) {
if (pass == RenderPass::Transparent) {
r->render(info);
}
} else if (pass == RenderPass::Opaque) {
r->render(info);
}
break;
}
}
}
/*========================================================================*/
CoordSet::CoordSet(PyMOLGlobals* G)
: CObjectState(G)
{
}
/*========================================================================*/
CoordSet::CoordSet(const CoordSet& cs)
: CObjectState(cs)
{
this->Obj = cs.Obj;
this->Coord = cs.Coord;
this->IdxToAtm = cs.IdxToAtm;
this->NIndex = cs.NIndex;
std::copy(std::begin(cs.Rep), std::end(cs.Rep), std::begin(this->Rep));
std::copy(std::begin(cs.Active), std::end(cs.Active), std::begin(this->Active));
this->NTmpBond = cs.NTmpBond;
this->NTmpLinkBond = cs.NTmpLinkBond;
/* deep copy & return ptr to new symmetry */
if(cs.Symmetry) {
this->Symmetry = std::make_unique<CSymmetry>(*cs.Symmetry);
}
std::copy(std::begin(cs.Name), std::end(cs.Name), std::begin(this->Name));
this->PeriodicBoxType = cs.PeriodicBoxType;
this->tmp_index = cs.tmp_index;
this->Coord2IdxReq = cs.Coord2IdxReq;
this->Coord2IdxDiv = cs.Coord2IdxDiv;
this->objMolOpInvalidated = cs.objMolOpInvalidated;
// copy VLAs
this->RefPos = cs.RefPos;
this->AtmToIdx = cs.AtmToIdx;
UtilZeroMem(this->Rep, sizeof(::Rep *) * cRepCnt);
if (cs.prop_id) {
PropertyCheckUniqueID(G, this);
PropertyCopyProperties(G, cs.prop_id, this->prop_id);
}
}
/*========================================================================*/
/**
* Sets the number of atoms with coordinates.
*
* @post NIndex updated
* @post All relevant arrays resized
*/
void CoordSet::setNIndex(unsigned nindex)
{
NIndex = nindex;
IdxToAtm.resize(nindex);
if (nindex == 0) {
return;
}
auto const idx = nindex - 1;
Coord.reserve(nindex * 3);
if (has_any_atom_state_settings()) {
atom_state_setting_id.check(idx);
}
if (RefPos) {
RefPos.check(idx);
}
}
/**
* Sets the number of atoms.
*
* For non-discrete objects:
* - Extends AtmToIdx (appends -1s)
*
* For discrete objects:
* - Converts this coordset to discrete if necessary
* - Calls ObjectMolecule::setNDiscrete()
*
* @pre NIndex and IdxToAtm are valid
*/
int CoordSet::extendIndices(int nAtom)
{
bool ok = true;
if (Obj->DiscreteFlag) {
ok = Obj->setNDiscrete(nAtom);
// convert to discrete if necessary
if (!AtmToIdx.empty()) {
AtmToIdx.clear();
if (ok) {
for (int a = 0; a < NIndex; a++) {
auto const b = IdxToAtm[a];
Obj->DiscreteAtmToIdx[b] = a;
Obj->DiscreteCSet[b] = this;
}
}
}
} else {
const auto NAtIndex = AtmToIdx.size();
assert(NAtIndex <= nAtom);
if (NAtIndex < nAtom) {
AtmToIdx.resize(nAtom);
if (ok && nAtom) {
for (int a = NAtIndex; a < nAtom; a++)
AtmToIdx[a] = -1;
}
}
}
return ok;
}
/*========================================================================*/
/**
* Set up for simple case where 1 = 1, etc.
*/
void CoordSet::enumIndices()
{
AtmToIdx.resize(NIndex);
IdxToAtm.resize(NIndex);
if (NIndex) {
for (int a = 0; a < NIndex; ++a) {
AtmToIdx[a] = a;
IdxToAtm[a] = a;
}
}
}
/*========================================================================*/
CoordSet::~CoordSet()
{
if (prop_id) {
PropertyUniqueDetachChain(G, prop_id);
prop_id = 0;
}
if (has_any_atom_state_settings()) {
for (int a = 0; a < NIndex; ++a) {
if (has_atom_state_settings(a)) {
SettingUniqueDetachChain(G, atom_state_setting_id[a]);
}
}
}
for (int a = 0; a < cRepCnt; ++a) {
delete Rep[a];
}
MapFree(Coord2Idx);
CGOFree(SculptCGO);
CGOFree(SculptShaderCGO);
}
void LabPosTypeCopy(const LabPosType * src, LabPosType * dst){
dst->mode = src->mode;
copy3f(src->pos, dst->pos);
copy3f(src->offset, dst->offset);
}
void RefPosTypeCopy(const RefPosType * src, RefPosType * dst){
copy3f(src->coord, dst->coord);
dst->specified = src->specified;
}
#ifndef _PYMOL_NOPY
int CoordSetSetSettingFromPyObject(PyMOLGlobals * G, CoordSet *cs, int at, int setting_id, PyObject *val){
if (val == Py_None)
val = nullptr;
if (!val && !cs->has_atom_state_settings(at)) {
return true;
}
CoordSetCheckUniqueID(G, cs, at);
return SettingUniqueSetPyObject(G, cs->atom_state_setting_id[at], setting_id, val);
}
#endif
int CoordSetCheckSetting(PyMOLGlobals * G, CoordSet *cs, int at, int setting_id){
if (!cs->has_atom_state_settings(at)) {
return 0;
}
return SettingUniqueCheck(G, cs->atom_state_setting_id[at], setting_id);
}
PyObject *SettingGetIfDefinedPyObject(PyMOLGlobals * G, CoordSet *cs, int at, int setting_id){
if (cs->has_atom_state_settings(at)) {
return SettingUniqueGetPyObject(G, cs->atom_state_setting_id[at], setting_id);
}
return nullptr;
}
int CoordSetCheckUniqueID(PyMOLGlobals * G, CoordSet *I, int at){
if (!I->atom_state_setting_id){
I->atom_state_setting_id = pymol::vla<int>(I->NIndex);
}
if (!I->atom_state_setting_id[at]){
I->atom_state_setting_id[at] = AtomInfoGetNewUniqueID(G);
}
return I->atom_state_setting_id[at];
}
template <typename V>
void AtomStateGetSetting(ATOMSTATEGETSETTINGARGS, V * out) {
if (cs->has_atom_state_settings(idx) &&
SettingUniqueGetIfDefined(G, cs->atom_state_setting_id[idx], setting_id, out))
return;
if (AtomSettingGetIfDefined(G, ai, setting_id, out))
return;
*out = SettingGet<V>(*cs, setting_id);
}
template void AtomStateGetSetting(ATOMSTATEGETSETTINGARGS, int * out);
template void AtomStateGetSetting(ATOMSTATEGETSETTINGARGS, float * out);
template void AtomStateGetSetting(ATOMSTATEGETSETTINGARGS, const float ** out);
Reference in New Issue View Git Blame Copy Permalink