1
2 /*
3 A* -------------------------------------------------------------------
4 B* This file contains source code for the PyMOL computer program
5 C* copyright 1998-2000 by Warren Lyford Delano of DeLano Scientific.
6 D* -------------------------------------------------------------------
7 E* It is unlawful to modify or remove this copyright notice.
8 F* -------------------------------------------------------------------
9 G* Please see the accompanying LICENSE file for further information.
10 H* -------------------------------------------------------------------
11 I* Additional authors of this source file include:
12 -*
13 -*
14 -*
15 Z* -------------------------------------------------------------------
16 */
17 #include"os_python.h"
18 #include"os_numpy.h"
19
20 #include"os_std.h"
21
22 #include <algorithm>
23
24 #include"Base.h"
25 #include"OOMac.h"
26 #include"MemoryDebug.h"
27 #include"Err.h"
28 #include"Scene.h"
29 #include"CoordSet.h"
30 #include"Color.h"
31 #include"PConv.h"
32 #include"P.h"
33 #include"Matrix.h"
34 #include"Sphere.h"
35 #include"Util.h"
36 #include"Feedback.h"
37 #include"RepWireBond.h"
38 #include"RepCylBond.h"
39 #include"RepDot.h"
40 #include"RepMesh.h"
41 #include"RepSphere.h"
42 #include"RepSphereImmediate.h"
43 #include"RepRibbon.h"
44 #include"RepCartoon.h"
45 #include"RepSurface.h"
46 #include"RepLabel.h"
47 #include"RepNonbonded.h"
48 #include"RepNonbondedSphere.h"
49 #include"RepEllipsoid.h"
50
51 #include"PyMOLGlobals.h"
52 #include"PyMOLObject.h"
53 #include "Executive.h"
54 #include "Lex.h"
55
56 #ifdef _PYMOL_IP_PROPERTIES
57 #include "Property.h"
58 #endif
59
60
61 /*========================================================================*/
62 /*
63 * Get coordinate index for given atom index
64 */
atmToIdx(int atm) const65 int CoordSet::atmToIdx(int atm) const {
66 if (Obj->DiscreteFlag) {
67 if (this == Obj->DiscreteCSet[atm])
68 return Obj->DiscreteAtmToIdx[atm];
69 return -1;
70 }
71 return AtmToIdx[atm];
72 }
73
74 /*========================================================================*/
75 static char sATOM[] = "ATOM ";
76 static char sHETATM[] = "HETATM";
77
CoordSetValidateRefPos(CoordSet * I)78 int CoordSetValidateRefPos(CoordSet * I)
79 {
80 if(I->RefPos) {
81 VLACheck(I->RefPos, RefPosType, I->NIndex);
82 return true;
83 } else {
84 int ok = true && (I->RefPos = pymol::vla<RefPosType>(I->NIndex));
85 if(ok) {
86 int a;
87 for(a = 0; a < I->NIndex; a++) {
88 const float* src = I->coordPtr(a);
89 copy3f(src, I->RefPos[a].coord);
90 I->RefPos[a].specified = true;
91 }
92 }
93 return ok;
94 }
95 }
96
97
98 /*========================================================================*/
BondCompare(BondType * a,BondType * b)99 int BondCompare(BondType * a, BondType * b)
100 {
101 int ai0 = a->index[0];
102 int bi0 = b->index[0];
103 if(ai0 == bi0) {
104 int ai1 = a->index[1];
105 int bi1 = b->index[1];
106 if(ai1 == bi1) {
107 return 0;
108 } else if(ai1 > bi1) {
109 return 1;
110 } else {
111 return -1;
112 }
113 } else if(ai0 > bi0) {
114 return 1;
115 } else {
116 return -1;
117 }
118 }
119
120
121 /*========================================================================*/
BondInOrder(BondType * a,int b1,int b2)122 int BondInOrder(BondType * a, int b1, int b2)
123 {
124 return (BondCompare(a + b1, a + b2) <= 0);
125 }
126
CoordSetFromPyList(PyMOLGlobals * G,PyObject * list,CoordSet ** cs)127 int CoordSetFromPyList(PyMOLGlobals * G, PyObject * list, CoordSet ** cs)
128 {
129 CoordSet *I = NULL;
130 int ok = true;
131 int ll = 0;
132
133 if(*cs) {
134 (*cs)->fFree();
135 *cs = NULL;
136 }
137
138 if(list == Py_None) { /* allow None for CSet */
139 *cs = NULL;
140 } else {
141
142 if(ok)
143 I = CoordSetNew(G);
144 if(ok)
145 ok = (I != NULL);
146 if(ok)
147 ok = (list != NULL);
148 if(ok)
149 ok = PyList_Check(list);
150 if(ok)
151 ll = PyList_Size(list);
152 /* TO SUPPORT BACKWARDS COMPATIBILITY...
153 Always check ll when adding new PyList_GetItem's */
154 if(ok)
155 ok = PConvPyIntToInt(PyList_GetItem(list, 0), &I->NIndex);
156 if(ok)
157 ok = PConvPyIntToInt(PyList_GetItem(list, 1), &I->NAtIndex);
158 if(ok)
159 ok = PConvPyListToFloatVLA(PyList_GetItem(list, 2), &I->Coord);
160 if(ok)
161 ok = PConvPyListToIntVLA(PyList_GetItem(list, 3), &I->IdxToAtm);
162 if(ok && (ll > 5))
163 ok = CPythonVal_PConvPyStrToStr_From_List(G, list, 5, I->Name, sizeof(WordType));
164 if(ok && (ll > 6)){
165 CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 6);
166 ok = ObjectStateFromPyList(G, val, I);
167 CPythonVal_Free(val);
168 }
169 if(ok && (ll > 7)){
170 CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 7);
171 I->Setting = SettingNewFromPyList(G, val);
172 CPythonVal_Free(val);
173 }
174 if(ok && (ll > 8)){
175 CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 8);
176 ok = CPythonVal_PConvPyListToLabPosVLA(G, val, &I->LabPos);
177 CPythonVal_Free(val);
178 }
179
180 #ifdef _PYMOL_IP_PROPERTIES
181 #endif
182
183 if(ok && (ll > 10)){
184 CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 10);
185 if (!CPythonVal_IsNone(val))
186 I->SculptCGO = CGONewFromPyList(G, val, 0, 1);
187 else {
188 I->SculptShaderCGO = I->SculptCGO = NULL;
189 }
190 CPythonVal_Free(val);
191 }
192 if(ok){
193 if(ll > 11){
194 // load in atom-state settings
195 CPythonVal *val = CPythonVal_PyList_GetItem(G, list, 11);
196 if (!CPythonVal_IsNone(val)){
197 int a;
198 I->atom_state_setting_id = pymol::vla<int>(I->NIndex);
199 I->has_atom_state_settings = pymol::vla<char>(I->NIndex);
200 for (a=0; a<I->NIndex; a++){
201 CPythonVal *val2 = CPythonVal_PyList_GetItem(G, val, a);
202 if (!CPythonVal_IsNone(val2)){
203 CPythonVal_PConvPyIntToInt(val2, &I->atom_state_setting_id[a]);
204 I->has_atom_state_settings[a] = (I->atom_state_setting_id[a]) ? 1 : 0;
205 if (I->atom_state_setting_id[a]) {
206 I->atom_state_setting_id[a] = SettingUniqueConvertOldSessionID(G, I->atom_state_setting_id[a]);
207 }
208 }
209 CPythonVal_Free(val2);
210 }
211 } else {
212 I->atom_state_setting_id = NULL;
213 I->has_atom_state_settings = NULL;
214 }
215 CPythonVal_Free(val);
216 } else {
217 // check I->LabPos.offset to see if its set
218 if (I->LabPos){
219 int a;
220 for (a=0; a<I->NIndex; a++){
221 if(length3f(I->LabPos[a].offset) > R_SMALL4) {
222 SettingSet(cSetting_label_placement_offset, I->LabPos[a].offset, I, a);
223 }
224 }
225 }
226 }
227 }
228 if(!ok) {
229 if(I)
230 I->fFree();
231 *cs = NULL;
232 } else {
233 *cs = I;
234 }
235 }
236 return (ok);
237 }
238
239 /*
240 * Coord set as numpy array
241 */
CoordSetAsNumPyArray(CoordSet * cs,short copy)242 PyObject *CoordSetAsNumPyArray(CoordSet * cs, short copy)
243 {
244 #ifndef _PYMOL_NUMPY
245 PRINTFB(cs->G, FB_CoordSet, FB_Errors)
246 "No numpy support\n" ENDFB(cs->G);
247 return NULL;
248 #else
249
250 PyObject *result = NULL;
251 const int base_size = sizeof(float);
252 int typenum = -1;
253 npy_intp dims[2] = {0, 3};
254
255 import_array1(NULL);
256
257 switch(base_size) {
258 case 4: typenum = NPY_FLOAT32; break;
259 case 8: typenum = NPY_FLOAT64; break;
260 }
261
262 if(typenum == -1) {
263 printf("error: no typenum for float size %d\n", base_size);
264 return NULL;
265 }
266
267 dims[0] = cs->NIndex;
268
269 if(copy) {
270 if((result = PyArray_SimpleNew(2, dims, typenum)))
271 memcpy(PyArray_DATA((PyArrayObject *)result), cs->Coord, cs->NIndex * 3 * base_size);
272 } else {
273 result = PyArray_SimpleNewFromData(2, dims, typenum, cs->Coord.data());
274 }
275
276 return result;
277 #endif
278 }
279
CoordSetAsPyList(CoordSet * I)280 PyObject *CoordSetAsPyList(CoordSet * I)
281 {
282 PyObject *result = NULL;
283
284 if(I) {
285 auto G = I->G;
286 int pse_export_version = SettingGet<float>(G, cSetting_pse_export_version) * 1000;
287 bool dump_binary = SettingGet<bool>(G, cSetting_pse_binary_dump) && (!pse_export_version || pse_export_version >= 1765);
288 result = PyList_New(12);
289 PyList_SetItem(result, 0, PyInt_FromLong(I->NIndex));
290 PyList_SetItem(result, 1, PyInt_FromLong(I->NAtIndex));
291 PyList_SetItem(result, 2, PConvFloatArrayToPyList(I->Coord, I->NIndex * 3, dump_binary));
292 PyList_SetItem(result, 3, PConvIntArrayToPyList(I->IdxToAtm, I->NIndex, dump_binary));
293 if(I->AtmToIdx
294 && pse_export_version < 1770)
295 PyList_SetItem(result, 4, PConvIntArrayToPyList(I->AtmToIdx, I->NAtIndex, dump_binary));
296 else
297 PyList_SetItem(result, 4, PConvAutoNone(NULL));
298 PyList_SetItem(result, 5, PyString_FromString(I->Name));
299 PyList_SetItem(result, 6, ObjectStateAsPyList(I));
300 PyList_SetItem(result, 7, SettingAsPyList(I->Setting));
301 PyList_SetItem(result, 8, PConvLabPosVLAToPyList(I->LabPos, I->NIndex));
302
303 PyList_SetItem(result, 9,
304 #ifdef _PYMOL_IP_PROPERTIES
305 #endif
306 PConvAutoNone(Py_None));
307
308 if(I->SculptCGO) {
309 PyList_SetItem(result, 10, CGOAsPyList(I->SculptCGO));
310 } else {
311 PyList_SetItem(result, 10, PConvAutoNone(NULL));
312 }
313 if (I->has_atom_state_settings){
314 int a;
315 PyObject *settings_list = NULL;
316 settings_list = PyList_New(I->NIndex);
317 for (a=0; a<I->NIndex; a++){
318 if (I->has_atom_state_settings[a]){
319 PyList_SetItem(settings_list, a, PyInt_FromLong(I->atom_state_setting_id[a]));
320 } else {
321 PyList_SetItem(settings_list, a, PConvAutoNone(NULL));
322 }
323 }
324 PyList_SetItem(result, 11, settings_list);
325 } else {
326 PyList_SetItem(result, 11, PConvAutoNone(NULL));
327 }
328 /* TODO symmetry, spheroid, periodic box ... */
329 }
330 return (PConvAutoNone(result));
331 }
332
CoordSetAdjustAtmIdx(CoordSet * I,int * lookup,int nAtom)333 void CoordSetAdjustAtmIdx(CoordSet * I, int *lookup, int nAtom)
334
335 /* performs second half of removal */
336 {
337 /* NOTE: only works in a compressive mode, where lookup[a]<=a */
338 int a, a0, atm;
339 char *new_has_atom_state_settings_by_atom = NULL;
340 int *new_atom_state_setting_id_by_atom = NULL;
341 int nIndex = I->NIndex;
342
343 if (I->has_atom_state_settings){
344 new_has_atom_state_settings_by_atom = VLACalloc(char, nIndex);
345 new_atom_state_setting_id_by_atom = VLACalloc(int, nIndex);
346 }
347
348 for(a = 0; a < nIndex; a++) {
349 atm = I->IdxToAtm[a];
350 a0 = lookup[atm];
351 if (a0 < 0){
352 if (I->has_atom_state_settings){
353 if (I->has_atom_state_settings[a]){
354 SettingUniqueDetachChain(I->G, I->atom_state_setting_id[a]);
355 I->has_atom_state_settings[a] = 0;
356 I->atom_state_setting_id[a] = 0;
357 }
358 }
359 } else if (new_has_atom_state_settings_by_atom){
360 new_has_atom_state_settings_by_atom[a0] = I->has_atom_state_settings[a];
361 new_atom_state_setting_id_by_atom[a0] = I->atom_state_setting_id[a];
362 }
363 }
364
365 if (I->AtmToIdx){
366 for(a = 0; a < I->NAtIndex; a++) {
367 a0 = lookup[a];
368 if(a0 >= 0) {
369 I->AtmToIdx[a0] = I->AtmToIdx[a];
370 }
371 }
372 }
373 I->NAtIndex = nAtom;
374 if (I->AtmToIdx){
375 VLASize(I->AtmToIdx, int, nAtom);
376 }
377 for(a = 0; a < I->NIndex; a++) {
378 atm = I->IdxToAtm[a] = lookup[I->IdxToAtm[a]];
379 if (new_has_atom_state_settings_by_atom){
380 I->has_atom_state_settings[a] = new_has_atom_state_settings_by_atom[atm];
381 I->atom_state_setting_id[a] = new_atom_state_setting_id_by_atom[atm];
382 }
383 }
384 if (new_has_atom_state_settings_by_atom){
385 VLAFreeP(new_has_atom_state_settings_by_atom);
386 VLAFreeP(new_atom_state_setting_id_by_atom);
387 }
388 PRINTFD(I->G, FB_CoordSet)
389 " CoordSetAdjustAtmIdx-Debug: leaving... NAtIndex: %d NIndex %d\n",
390 I->NAtIndex, I->NIndex ENDFD;
391
392 }
393
CoordSetCopy(const CoordSet * src)394 CoordSet* CoordSetCopy(const CoordSet* src)
395 {
396 if(!src) {
397 return nullptr;
398 }
399 return new CoordSet(*src);
400 }
401
402
403 /*========================================================================*/
CoordSetMerge(ObjectMolecule * OM,CoordSet * I,CoordSet * cs)404 int CoordSetMerge(ObjectMolecule *OM, CoordSet * I, CoordSet * cs)
405 { /* must be non-overlapping */
406 int nIndex;
407 int a, i0;
408 int ok = true;
409 /* calculate new size and make room for new data */
410 nIndex = I->NIndex + cs->NIndex;
411 VLASize(I->IdxToAtm, int, nIndex);
412 CHECKOK(ok, I->IdxToAtm);
413 if (ok)
414 VLACheck(I->Coord, float, nIndex * 3);
415 CHECKOK(ok, I->Coord);
416 if (ok){
417 for(a = 0; a < cs->NIndex; a++) {
418 i0 = a + I->NIndex;
419 I->IdxToAtm[i0] = cs->IdxToAtm[a];
420 if (OM->DiscreteFlag){
421 int idx = cs->IdxToAtm[a];
422 OM->DiscreteAtmToIdx[idx] = i0;
423 OM->DiscreteCSet[idx] = I;
424 } else {
425 I->AtmToIdx[cs->IdxToAtm[a]] = i0;
426 }
427 copy3f(cs->coordPtr(a), I->coordPtr(i0));
428 }
429 }
430 if (ok){
431 if(cs->LabPos) {
432 if(!I->LabPos)
433 I->LabPos = pymol::vla<LabPosType>(nIndex);
434 else
435 VLACheck(I->LabPos, LabPosType, nIndex);
436 if(I->LabPos) {
437 UtilCopyMem(I->LabPos + I->NIndex, cs->LabPos, sizeof(LabPosType) * cs->NIndex);
438 }
439 } else if(I->LabPos) {
440 VLACheck(I->LabPos, LabPosType, nIndex);
441 }
442 }
443 if (ok){
444 if(cs->RefPos) {
445 if(!I->RefPos)
446 I->RefPos = pymol::vla<RefPosType>(nIndex);
447 else
448 VLACheck(I->RefPos, RefPosType, nIndex);
449 if(I->RefPos) {
450 UtilCopyMem(I->RefPos + I->NIndex, cs->RefPos, sizeof(RefPosType) * cs->NIndex);
451 }
452 } else if(I->RefPos) {
453 VLACheck(I->RefPos, RefPosType, nIndex);
454 }
455 I->invalidateRep(cRepAll, cRepInvAll);
456 }
457 I->NIndex = nIndex;
458
459 return ok;
460 }
461
462
463 /*========================================================================*/
CoordSetPurge(CoordSet * I)464 void CoordSetPurge(CoordSet * I)
465
466 /* performs first half of removal */
467 {
468 auto G = I->G;
469 int offset = 0;
470 int a, a1, ao;
471 AtomInfoType *ai;
472 ObjectMolecule *obj;
473 float *c0, *c1;
474 LabPosType *l0, *l1;
475 RefPosType *r0, *r1;
476 int *atom_state0, *atom_state1;
477 char *has_atom_state0, *has_atom_state1;
478 obj = I->Obj;
479
480 PRINTFD(G, FB_CoordSet)
481 " CoordSetPurge-Debug: entering..." ENDFD;
482
483 c0 = c1 = I->Coord.data();
484 r0 = r1 = I->RefPos.data();
485 l0 = l1 = I->LabPos.data();
486 atom_state0 = atom_state1 = I->atom_state_setting_id.data();
487 has_atom_state0 = has_atom_state1 = I->has_atom_state_settings.data();
488
489 /* This loop slides down the atoms that are not deleted (deleteFlag)
490 it moves the Coord, RefPos, and LabPos */
491 for(a = 0; a < I->NIndex; a++) {
492 a1 = I->IdxToAtm[a];
493 ai = obj->AtomInfo + a1;
494 if(ai->deleteFlag) {
495 offset--;
496 c0 += 3;
497 if(l0)
498 l0++;
499 if(r0)
500 r0++;
501 if (has_atom_state0){
502 atom_state0++;
503 has_atom_state0++;
504 }
505 } else if(offset) {
506 ao = a + offset;
507 *(c1++) = *(c0++);
508 *(c1++) = *(c0++);
509 *(c1++) = *(c0++);
510 if(r1) {
511 *(r1++) = *(r0++);
512 }
513 if(l0) {
514 *(l1++) = *(l0++);
515 }
516 if (has_atom_state0){
517 *(atom_state1++) = *(atom_state0++);
518 *(has_atom_state1++) = *(has_atom_state0++);
519 }
520 if (I->AtmToIdx)
521 I->AtmToIdx[a1] = ao;
522 I->IdxToAtm[ao] = a1; /* no adjustment of these indexes yet... */
523 if (I->Obj->DiscreteFlag){
524 I->Obj->DiscreteAtmToIdx[a1] = ao;
525 I->Obj->DiscreteCSet[a1] = I;
526 }
527 } else {
528 c0 += 3;
529 c1 += 3;
530 if(r1) {
531 r0++;
532 r1++;
533 }
534 if(l0) {
535 l0++;
536 l1++;
537 }
538 if (has_atom_state0){
539 atom_state0++; atom_state1++;
540 has_atom_state0++; has_atom_state1++;
541 }
542 }
543 }
544 if(offset) {
545 /* If there were deleted atoms, (offset < 0), then
546 re-adjust the array sizes */
547 I->NIndex += offset;
548 VLASize(I->Coord, float, I->NIndex * 3);
549 if(I->LabPos) {
550 VLASize(I->LabPos, LabPosType, I->NIndex);
551 }
552 if(I->RefPos) {
553 VLASize(I->RefPos, RefPosType, I->NIndex);
554 }
555 if(I->has_atom_state_settings) {
556 VLASize(I->has_atom_state_settings, char, I->NIndex);
557 VLASize(I->atom_state_setting_id, int, I->NIndex);
558 }
559 VLASize(I->IdxToAtm, int, I->NIndex);
560 PRINTFD(G, FB_CoordSet)
561 " CoordSetPurge-Debug: I->IdxToAtm shrunk to %d\n", I->NIndex ENDFD;
562 I->invalidateRep(cRepAll, cRepInvAtoms); /* this will free Color */
563 }
564 PRINTFD(G, FB_CoordSet)
565 " CoordSetPurge-Debug: leaving NAtIndex %d NIndex %d...\n",
566 I->NAtIndex, I->NIndex ENDFD;
567
568 #ifdef _PYMOL_IP_PROPERTIES
569 #endif
570 }
571
572
573 /*========================================================================*/
CoordSetTransformAtomTTTf(CoordSet * I,int at,const float * TTT)574 int CoordSetTransformAtomTTTf(CoordSet * I, int at, const float *TTT)
575 {
576 int a1 = I->atmToIdx(at);
577 if(a1 < 0)
578 return false;
579
580 auto* v1 = I->coordPtr(a1);
581 MatrixTransformTTTfN3f(1, v1, TTT, v1);
582 return true;
583 }
584
585
586 /*========================================================================*/
CoordSetTransformAtomR44f(CoordSet * I,int at,const float * matrix)587 int CoordSetTransformAtomR44f(CoordSet * I, int at, const float *matrix)
588 {
589 int a1 = I->atmToIdx(at);
590 if(a1 < 0)
591 return false;
592
593 auto* v1 = I->coordPtr(a1);
594 MatrixTransformR44fN3f(1, v1, matrix, v1);
595 return true;
596 }
597
598
599 /*========================================================================*/
CoordSetRecordTxfApplied(CoordSet * I,const float * matrix,int homogenous)600 void CoordSetRecordTxfApplied(CoordSet * I, const float *matrix, int homogenous)
601 {
602 double temp[16];
603
604 if(!homogenous) {
605 convertTTTfR44d(matrix, temp);
606 } else {
607 convert44f44d(matrix, temp);
608 }
609
610 ObjectStateLeftCombineMatrixR44d(I, temp);
611 }
612
613
614 /*========================================================================*/
CoordSetMoveAtom(CoordSet * I,int at,const float * v,int mode)615 int CoordSetMoveAtom(CoordSet * I, int at, const float *v, int mode)
616 {
617 int a1 = I->atmToIdx(at);
618 if(a1 < 0)
619 return false;
620
621 auto* v1 = I->coordPtr(a1);
622 if(mode) {
623 add3f(v, v1, v1);
624 } else {
625 copy3f(v, v1);
626 }
627 return true;
628 }
629
630
631 /*========================================================================*/
CoordSetMoveAtomLabel(CoordSet * I,int at,const float * v,const float * diff)632 int CoordSetMoveAtomLabel(CoordSet * I, int at, const float *v, const float *diff)
633 {
634 auto G = I->G;
635 ObjectMolecule *obj = I->Obj;
636 int a1 = I->atmToIdx(at);
637 int result = 0;
638
639 /* if label is valid, get the label offset
640 * and set the new position relative to that */
641 if(a1 >= 0) {
642 float at_offset[3];
643 const float * at_offset_ptr;
644 int at_label_relative_mode = 0;
645 AtomInfoType *ai = obj->AtomInfo + at;
646
647 AtomStateGetSetting_i(G, obj, I, a1, ai, cSetting_label_relative_mode, &at_label_relative_mode);
648 switch (at_label_relative_mode){
649 case 0:
650 AtomStateGetSetting(G, obj, I, a1, ai, cSetting_label_placement_offset, &at_offset_ptr);
651 add3f(v, at_offset_ptr, at_offset);
652 SettingSet(cSetting_label_placement_offset, at_offset, I, a1);
653 break;
654 case 1: // screen relative
655 case 2: // screen pixel space
656 {
657 float voff[3];
658 int width, height;
659 SceneGetWidthHeight(G, &width, &height);
660 if (at_label_relative_mode==1){
661 voff[0] = 2.f * diff[0] / width;
662 voff[1] = 2.f * diff[1] / height;
663 } else {
664 voff[0] = diff[0];
665 voff[1] = diff[1];
666 }
667 voff[2] = 0.f;
668 AtomStateGetSetting(G, obj, I, a1, ai, cSetting_label_screen_point, &at_offset_ptr);
669 add3f(voff, at_offset_ptr, at_offset);
670 SettingSet(cSetting_label_screen_point, at_offset, I, a1);
671 }
672 break;
673 }
674 }
675
676 return (result);
677 }
678
679
680 /*========================================================================*/
CoordSetGetAtomVertex(const CoordSet * I,int at,float * v)681 int CoordSetGetAtomVertex(const CoordSet * I, int at, float *v)
682 {
683 int a1 = I->atmToIdx(at);
684
685 if(a1 < 0)
686 return false;
687
688 copy3f(I->coordPtr(a1), v);
689 return true;
690 }
691
692
693 /*========================================================================*/
CoordSetGetAtomTxfVertex(const CoordSet * I,int at,float * v)694 int CoordSetGetAtomTxfVertex(const CoordSet * I, int at, float *v)
695 {
696 ObjectMolecule *obj = I->Obj;
697 int a1 = I->atmToIdx(at);
698
699 if(a1 < 0)
700 return false;
701
702 copy3f(I->coordPtr(a1), v);
703
704 /* apply state transformation */
705 if (!I->Matrix.empty() && SettingGet<int>(I->G, obj->Setting, I->Setting,
706 cSetting_matrix_mode) > 0) {
707 transform44d3f(I->Matrix.data(), v, v);
708 }
709
710 /* object transformation */
711 if(obj->TTTFlag) {
712 transformTTT44f3f(obj->TTT, v, v);
713 }
714
715 return true;
716 }
717
718
719 /*========================================================================*/
CoordSetSetAtomVertex(CoordSet * I,int at,const float * v)720 int CoordSetSetAtomVertex(CoordSet * I, int at, const float *v)
721 {
722 int a1 = I->atmToIdx(at);
723
724 if(a1 < 0)
725 return false;
726
727 copy3f(v, I->coordPtr(a1));
728 return true;
729 }
730
731
732 /*========================================================================*/
CoordSetRealToFrac(CoordSet * I,const CCrystal * cryst)733 void CoordSetRealToFrac(CoordSet * I, const CCrystal * cryst)
734 {
735 int a;
736 float* v = I->Coord.data();
737 for(a = 0; a < I->NIndex; a++) {
738 transform33f3f(cryst->RealToFrac, v, v);
739 v += 3;
740 }
741 }
742
743
744 /*========================================================================*/
CoordSetTransform44f(CoordSet * I,const float * mat)745 void CoordSetTransform44f(CoordSet * I, const float *mat)
746 {
747 int a;
748 float* v = I->Coord.data();
749 for(a = 0; a < I->NIndex; a++) {
750 transform44f3f(mat, v, v);
751 v += 3;
752 }
753 }
754
755
756 /*========================================================================*/
757
CoordSetTransform33f(CoordSet * I,const float * mat)758 void CoordSetTransform33f(CoordSet * I, const float *mat)
759 {
760 int a;
761 float* v = I->Coord.data();
762 for(a = 0; a < I->NIndex; a++) {
763 transform33f3f(mat, v, v);
764 v += 3;
765 }
766 }
767
768
769 /*========================================================================*/
CoordSetGetAverage(const CoordSet * I,float * v0)770 void CoordSetGetAverage(const CoordSet * I, float *v0)
771 {
772 int a;
773 double accum[3];
774 if(I->NIndex) {
775 const float* v = I->Coord.data();
776 accum[0] = *(v++);
777 accum[1] = *(v++);
778 accum[2] = *(v++);
779 for(a = 1; a < I->NIndex; a++) {
780 accum[0] += *(v++);
781 accum[1] += *(v++);
782 accum[2] += *(v++);
783 }
784 v0[0] = (float) (accum[0] / I->NIndex);
785 v0[1] = (float) (accum[1] / I->NIndex);
786 v0[2] = (float) (accum[2] / I->NIndex);
787 }
788 }
789
790
791 /*========================================================================*/
CoordSetFracToReal(CoordSet * I,const CCrystal * cryst)792 void CoordSetFracToReal(CoordSet * I, const CCrystal * cryst)
793 {
794 int a;
795 float* v = I->Coord.data();
796 for(a = 0; a < I->NIndex; a++) {
797 transform33f3f(cryst->FracToReal, v, v);
798 v += 3;
799 }
800 }
801
802 /*
803 * Apply the `sca` (SCALEn) transformation to transform to fractional space,
804 * and then the crystals `FracToReal` transformation to transform back to
805 * cartesian space.
806 *
807 * Don't do anything if pdb_insure_orthogonal=off.
808 *
809 * Don't do anything if SCALEn or CRYST1 look bogus. There is a number of
810 * structures in the PDB which have meaningless values for those.
811 *
812 * Without this, creating symmetry mates might produce wrong results.
813 */
CoordSetInsureOrthogonal(PyMOLGlobals * G,CoordSet * cset,const float * sca,const CCrystal * cryst,bool quiet)814 bool CoordSetInsureOrthogonal(PyMOLGlobals * G,
815 CoordSet * cset, // coord set to modify
816 const float * sca, // 4x4 SCALE
817 const CCrystal *cryst,
818 bool quiet)
819 {
820 if (!SettingGetGlobal_b(G, cSetting_pdb_insure_orthogonal))
821 return false;
822
823 if (!cryst)
824 cryst = &cset->Symmetry->Crystal;
825
826 const float * r2f = cryst->RealToFrac;
827
828 // are the matrices sufficiently close to be the same?
829 if (!sca[3] && !sca[7] && !sca[11] &&
830 is_allclosef(3, r2f, 3, sca, 4, R_SMALL4)) {
831 return false;
832 }
833
834 // is the cell a orthogonal 1x1x1? If so, then it should probably be ignored...
835 // is SCALEn the identity matrix? If so, then it should probably be ignored...
836 if (is_identityf(3, r2f, R_SMALL4) ||
837 is_identityf(4, sca, R_SMALL4)) {
838 PRINTFB(G, FB_ObjectMolecule, FB_Blather)
839 " ObjectMolReadPDBStr: ignoring SCALEn (identity matrix).\n" ENDFB(G);
840 return false;
841 }
842
843 // is SCALEn invalid? If so, then it should definitely be ignored...
844 if (determinant33f(sca, 4) < R_SMALL8 ||
845 determinant33f(r2f, 3) < R_SMALL8) {
846 PRINTFB(G, FB_ObjectMolecule, FB_Blather)
847 " ObjectMolReadPDBStr: ignoring SCALEn (invalid matrix).\n" ENDFB(G);
848 return false;
849 }
850
851 PRINTFB(G, FB_ObjectMolecule, quiet ? FB_Blather : FB_Actions)
852 " ObjectMolecule: using SCALEn to compute orthogonal coordinates.\n"
853 ENDFB(G);
854
855 CoordSetTransform44f(cset, sca);
856 CoordSetFracToReal(cset, cryst);
857
858 return true;
859 }
860
861 /*========================================================================*/
862 /* Rotates the ANISOU vector
863 *
864 * matrix: flat 4x4, but only rotation (upper left 3x3) is considered
865 * anisou: has 6 elements (of symmetric 3x3) and will be rotated in-place
866 */
RotateU(const double * matrix,float * anisou)867 bool RotateU(const double *matrix, float *anisou)
868 {
869 int i, j, k;
870 float Re[3][3];
871 double e_val[3], e_vec[3][3];
872 double U[3][3] = {
873 { anisou[0], anisou[3], anisou[4] },
874 { anisou[3], anisou[1], anisou[5] },
875 { anisou[4], anisou[5], anisou[2] },
876 };
877
878 // e_val, e_vec = linalg.eigh(U)
879 if(!xx_matrix_jacobi_solve(*e_vec, e_val, &i, *U, 3))
880 return false;
881
882 // Re = dot(matrix[:3,:3], e_vec)
883 for (i = 0; i < 3; i++)
884 for (j = 0; j < 3; j++) {
885 Re[i][j] = 0.0;
886 for (k = 0; k < 3; k++)
887 Re[i][j] += matrix[i * 4 + k] * e_vec[k][j];
888 }
889
890 // U = dot(Re * e_val, Re.T)
891 for (i = 0; i < 3; i++)
892 for (j = 0; j <= i; j++) {
893 U[i][j] = 0.0;
894 for (k = 0; k < 3; k++)
895 U[i][j] += Re[i][k] * e_val[k] * Re[j][k];
896 }
897
898 anisou[0] = U[0][0];
899 anisou[1] = U[1][1];
900 anisou[2] = U[2][2];
901 anisou[3] = U[1][0];
902 anisou[4] = U[2][0];
903 anisou[5] = U[2][1];
904
905 return true;
906 }
907
908 /*========================================================================*/
CoordSetAtomToPDBStrVLA(PyMOLGlobals * G,char ** charVLA,int * c,const AtomInfoType * ai,const float * v,int cnt,const PDBInfoRec * pdb_info,const double * matrix)909 void CoordSetAtomToPDBStrVLA(PyMOLGlobals * G, char **charVLA, int *c,
910 const AtomInfoType * ai,
911 const float *v, int cnt,
912 const PDBInfoRec * pdb_info,
913 const double *matrix)
914 /*
915 * v: 3x1 vertex in final output space
916 * matrix: 4x4 homogenous transformation matrix from model space to output
917 * space (view matrix * state matrix). Used for ANISOU.
918 */
919 {
920 char *aType;
921 AtomName name;
922 ResName resn;
923 lexidx_t chain;
924
925 char formalCharge[4];
926 int ignore_pdb_segi = SettingGetGlobal_b(G, cSetting_ignore_pdb_segi);
927 WordType x, y, z;
928
929 AtomInfoGetAlignedPDBResidueName(G, ai, resn);
930 AtomInfoGetAlignedPDBAtomName(G, ai, resn, name);
931
932 formalCharge[0] = 0;
933 if(SettingGetGlobal_b(G, cSetting_pdb_formal_charges)) {
934 if((ai->formalCharge > 0) && (ai->formalCharge < 10)) {
935 sprintf(formalCharge, "%d+", ai->formalCharge);
936 } else if((ai->formalCharge < 0) && (ai->formalCharge > -10)) {
937 sprintf(formalCharge, "%d-", -ai->formalCharge);
938 }
939 }
940
941 if(ai->hetatm)
942 aType = sHETATM;
943 else
944 aType = sATOM;
945
946 char inscode = ai->getInscode(true);
947
948 VLACheck(*charVLA, char, (*c) + 1000);
949
950 if(SettingGetGlobal_b(G, cSetting_pdb_retain_ids)) {
951 cnt = ai->id - 1;
952 }
953 if(cnt > 99998)
954 cnt = 99998;
955
956 if((!pdb_info) || (!pdb_info->is_pqr_file())) { /* relying upon short-circuit */
957 short linelen;
958 sprintf(x, "%8.3f", v[0]);
959 x[8] = 0;
960 sprintf(y, "%8.3f", v[1]);
961 y[8] = 0;
962 sprintf(z, "%8.3f", v[2]);
963 z[8] = 0;
964 linelen =
965 sprintf((*charVLA) + (*c),
966 "%6s%5i %-4s%1s%-4s%1.1s%4i%c %s%s%s%6.2f%6.2f %-4.4s%2s%2s\n", aType,
967 cnt + 1, name, ai->alt, resn, LexStr(G, ai->chain), ai->resv % 10000, inscode, x, y, z, ai->q, ai->b,
968 ignore_pdb_segi ? "" :
969 LexStr(G, ai->segi), ai->elem, formalCharge);
970 if(ai->anisou) {
971 // Columns 7 - 27 and 73 - 80 are identical to the corresponding ATOM/HETATM record.
972 char *atomline = (*charVLA) + (*c);
973 char *anisoline = atomline + linelen;
974 float anisou[6];
975 std::copy_n(ai->anisou, 6, anisou);
976
977 if(matrix && !RotateU(matrix, anisou)) {
978 PRINTFB(G, FB_CoordSet, FB_Errors) "RotateU failed\n" ENDFB(G);
979 return;
980 }
981 strncpy(anisoline + 6, atomline + 6, 22);
982 sprintf(anisoline + 28,
983 "%7.0f%7.0f%7.0f%7.0f%7.0f%7.0f",
984 anisou[0] * 1e4, anisou[1] * 1e4, anisou[2] * 1e4,
985 anisou[3] * 1e4, anisou[4] * 1e4, anisou[5] * 1e4);
986 strcpy(anisoline + 70, atomline + 70);
987 strncpy(anisoline, "ANISOU", 6);
988 (*c) += linelen;
989 }
990 (*c) += linelen;
991 } else {
992 Chain alt;
993 if(pdb_info->is_pqr_file() && pdb_info->pqr_workarounds) {
994 inscode = ' ';
995 chain = 0; /* no chain IDs */
996 alt[0] = 0; /* not alt conf identifiers */
997 } else {
998 alt[0] = ai->alt[0];
999 alt[1] = 0;
1000 chain = ai->chain;
1001 }
1002 sprintf(x, "%8.3f", v[0]);
1003 if(x[0] != 32)
1004 sprintf(x, " %7.2f", v[0]);
1005 x[8] = 0;
1006 sprintf(y, "%8.3f", v[1]);
1007 y[8] = 0;
1008 if(y[0] != 32)
1009 sprintf(y, " %7.2f", v[1]);
1010 y[8] = 0;
1011 sprintf(z, "%8.3f", v[2]);
1012 if(z[0] != 32)
1013 sprintf(z, " %7.2f", v[2]);
1014 z[8] = 0;
1015
1016 (*c) += sprintf((*charVLA) + (*c), "%6s%5i %-4s%1s%-4s%1.1s%4i%c %s%s%s %11.8f %7.3f\n",
1017 aType, cnt + 1, name, alt, resn,
1018 LexStr(G, chain), ai->resv, inscode, x, y, z, ai->partialCharge, ai->elec_radius);
1019 }
1020
1021 }
1022
1023
1024 /*========================================================================*/
CoordSetAtomToChemPyAtom(PyMOLGlobals * G,AtomInfoType * ai,const float * v,const float * ref,int index,const double * matrix)1025 PyObject *CoordSetAtomToChemPyAtom(PyMOLGlobals * G, AtomInfoType * ai, const float *v,
1026 const float *ref, int index, const double *matrix)
1027 {
1028 #ifdef _PYMOL_NOPY
1029 return NULL;
1030 #else
1031 PyObject *atom = PYOBJECT_CALLMETHOD(P_chempy, "Atom", "");
1032 if(!atom)
1033 ErrMessage(G, "CoordSetAtomToChemPyAtom", "can't create atom");
1034 else {
1035 float tmp_array[6] = { 0.f, 0.f, 0.f, 0.f, 0.f, 0.f };
1036
1037 if (ai->anisou) {
1038 std::copy_n(ai->anisou, 6, tmp_array);
1039 if (matrix)
1040 RotateU(matrix, tmp_array);
1041 }
1042
1043 PConvFloat3ToPyObjAttr(atom, "coord", v);
1044 if(ref)
1045 PConvFloat3ToPyObjAttr(atom, "ref_coord", ref);
1046 if (ai->name)
1047 PConvStringToPyObjAttr(atom, "name", LexStr(G, ai->name));
1048 PConvStringToPyObjAttr(atom, "symbol", ai->elem);
1049 // TODO defaults to UNK if (ai->resn)
1050 PConvStringToPyObjAttr(atom, "resn", LexStr(G, ai->resn));
1051 if (ai->inscode) {
1052 char ins_code[2] = {ai->inscode, 0};
1053 PConvStringToPyObjAttr(atom, "ins_code", ins_code);
1054 }
1055 if (ai->ssType[0])
1056 PConvStringToPyObjAttr(atom, "ss", ai->ssType);
1057 // TODO defaults to 1 if (ai->resv)
1058 PConvIntToPyObjAttr(atom, "resi_number", ai->resv);
1059 if (ai->stereo)
1060 PConvIntToPyObjAttr(atom, "stereo", ai->stereo);
1061 if (ai->chain)
1062 PConvStringToPyObjAttr(atom, "chain", LexStr(G, ai->chain));
1063 if(ai->alt[0])
1064 PConvStringToPyObjAttr(atom, "alt", ai->alt);
1065 if (ai->segi)
1066 PConvStringToPyObjAttr(atom, "segi", LexStr(G, ai->segi));
1067 if (ai->q != 1.)
1068 PConvFloatToPyObjAttr(atom, "q", ai->q);
1069 if (ai->b)
1070 PConvFloatToPyObjAttr(atom, "b", ai->b);
1071 if (ai->anisou)
1072 {
1073 {
1074 PyObject *tmp_obj = PConvFloatArrayToPyList(tmp_array, 6);
1075 if(tmp_obj) {
1076 PyObject_SetAttrString(atom, "u_aniso", tmp_obj);
1077 Py_XDECREF(tmp_obj);
1078 }
1079 }
1080 }
1081 PConvFloatToPyObjAttr(atom, "vdw", ai->vdw);
1082 if (ai->elec_radius)
1083 PConvFloatToPyObjAttr(atom, "elec_radius", ai->elec_radius);
1084 if (ai->partialCharge)
1085 PConvFloatToPyObjAttr(atom, "partial_charge", ai->partialCharge);
1086 if (ai->formalCharge)
1087 PConvIntToPyObjAttr(atom, "formal_charge", ai->formalCharge);
1088 // TODO customType=0 from most files
1089 if(ai->customType != -9999)
1090 PConvIntToPyObjAttr(atom, "numeric_type", ai->customType);
1091 if (ai->textType)
1092 PConvStringToPyObjAttr(atom, "text_type", LexStr(G, ai->textType));
1093 if (ai->custom)
1094 PConvStringToPyObjAttr(atom, "custom", LexStr(G, ai->custom));
1095
1096 PConvIntToPyObjAttr(atom, "hetatm", ai->hetatm);
1097 PConvIntToPyObjAttr(atom, "flags", ai->flags);
1098 PConvIntToPyObjAttr(atom, "id", ai->id); /* not necc. unique */
1099 PConvIntToPyObjAttr(atom, "index", index + 1); /* fragile */
1100
1101 #ifdef _PYMOL_IP_PROPERTIES
1102 #endif
1103
1104 }
1105 if(PyErr_Occurred())
1106 PyErr_Print();
1107 return (atom);
1108 #endif
1109 }
1110
1111
1112 /*========================================================================*/
invalidateRep(int type,int level)1113 void CoordSet::invalidateRep(int type, int level)
1114 {
1115 CoordSet * I = this;
1116 int a;
1117 if(level >= cRepInvVisib) {
1118 if (I->Obj)
1119 I->Obj->RepVisCacheValid = false;
1120 }
1121 /* graphical representations need redrawing */
1122 if(level == cRepInvVisib) {
1123 /* cartoon_side_chain_helper */
1124 if(SettingGet<bool>(G, I->Setting, I->Obj->Setting,
1125 cSetting_cartoon_side_chain_helper)) {
1126 if((type == cRepCyl) || (type == cRepLine) || (type == cRepSphere))
1127 invalidateRep(cRepCartoon, cRepInvVisib2);
1128 else if(type == cRepCartoon) {
1129 invalidateRep(cRepLine, cRepInvVisib2);
1130 invalidateRep(cRepCyl, cRepInvVisib2);
1131 invalidateRep(cRepSphere, cRepInvVisib2);
1132 }
1133 }
1134 /* ribbon_side_chain_helper */
1135 if(SettingGet<bool>(G, I->Setting, I->Obj->Setting,
1136 cSetting_ribbon_side_chain_helper)) {
1137 if((type == cRepCyl) || (type == cRepLine) || (type == cRepSphere))
1138 invalidateRep(cRepRibbon, cRepInvVisib2);
1139 else if(type == cRepRibbon) {
1140 invalidateRep(cRepLine, cRepInvVisib2);
1141 invalidateRep(cRepCyl, cRepInvVisib2);
1142 invalidateRep(cRepSphere, cRepInvVisib2);
1143 }
1144 }
1145 /* line_stick helper */
1146 if(SettingGet<bool>(G, I->Setting, I->Obj->Setting,
1147 cSetting_line_stick_helper)) {
1148 if(type == cRepCyl)
1149 invalidateRep(cRepLine, cRepInvVisib2);
1150 else if(type == cRepLine) {
1151 invalidateRep(cRepCyl, cRepInvVisib2);
1152 }
1153 }
1154 }
1155
1156 if(!I->Spheroid.empty())
1157 if (I->Spheroid.size() !=
1158 I->NAtIndex * GetSpheroidSphereRec(G)->nDot) {
1159 I->Spheroid.clear();
1160 I->SpheroidNormal.clear();
1161 }
1162
1163 /* invalidate basd on one representation, 'type' */
1164 for (RepIterator iter(G, type); iter.next(); ){
1165 int eff_level = level;
1166 a = iter.rep;
1167 if(level == cRepInvPick) {
1168 switch (a) {
1169 case cRepSurface:
1170 case cRepMesh:
1171 case cRepDot:
1172 /* skip the expensive to recompute, non-pickable
1173 representations */
1174 break;
1175 default: /* default behavior is to blow away the representation */
1176 eff_level = cRepInvRep;
1177 break;
1178 }
1179 }
1180 if(eff_level >= cRepInvVisib) /* make active if visibility has changed */
1181 I->Active[a] = true;
1182 if(I->Rep[a]) {
1183 if(I->Rep[a]->fInvalidate && (eff_level < cRepInvPurge))
1184 I->Rep[a]->fInvalidate(I->Rep[a], I, eff_level);
1185 else if(eff_level >= cRepInvExtColor) {
1186 I->Rep[a]->fFree(I->Rep[a]);
1187 I->Rep[a] = NULL;
1188 }
1189 }
1190 }
1191
1192 if(level >= cRepInvCoord) { /* if coordinates change, then this map becomes invalid */
1193 MapFree(I->Coord2Idx);
1194 I->Coord2Idx = NULL;
1195 ExecutiveInvalidateSelectionIndicatorsCGO(G);
1196 SceneInvalidatePicking(G);
1197 /* invalidate distances */
1198 }
1199
1200 #ifndef NO_MMLIBS
1201 if (level >= cRepInvProp) {
1202 if (validMMStereo == MMPYMOLX_PROP_STATE_AUTO)
1203 validMMStereo = MMPYMOLX_PROP_STATE_NULL;
1204 if (validTextType == MMPYMOLX_PROP_STATE_AUTO)
1205 validTextType = MMPYMOLX_PROP_STATE_NULL;
1206 }
1207 #endif
1208
1209 SceneChanged(G);
1210 }
1211
1212
1213 /*========================================================================*/
1214
1215 #define RepUpdateMacro(I,rep,new_fn,state) {\
1216 if(I->Active[rep]&&(!G->Interrupt)) {\
1217 if(!I->Rep[rep]) {\
1218 I->Rep[rep]=new_fn(I,state);\
1219 if(I->Rep[rep]){ \
1220 I->Rep[rep]->fNew=(struct Rep *(*)(struct CoordSet *,int state))new_fn;\
1221 SceneInvalidatePicking(G);\
1222 } else { \
1223 I->Active[rep] = false; \
1224 } \
1225 } else {\
1226 if(I->Rep[rep]->fUpdate)\
1227 I->Rep[rep] = I->Rep[rep]->fUpdate(I->Rep[rep],I,state,rep);\
1228 }\
1229 }\
1230 OrthoBusyFast(G,rep,cRepCnt);\
1231 }
1232
1233 /*========================================================================*/
update(int state)1234 void CoordSet::update(int state)
1235 {
1236 CoordSet * I = this;
1237 int a;
1238 assert(G == I->Obj->G);
1239
1240 PRINTFB(G, FB_CoordSet, FB_Blather) " CoordSetUpdate-Entered: object %s state %d cset %p\n",
1241 I->Obj->Name, state, (void *) I
1242 ENDFB(G);
1243
1244 OrthoBusyFast(G, 0, cRepCnt);
1245 RepUpdateMacro(I, cRepLine, RepWireBondNew, state);
1246 RepUpdateMacro(I, cRepCyl, RepCylBondNew, state);
1247 RepUpdateMacro(I, cRepDot, RepDotNew, state);
1248 RepUpdateMacro(I, cRepMesh, RepMeshNew, state);
1249 RepUpdateMacro(I, cRepSphere, RepSphereNew, state);
1250 RepUpdateMacro(I, cRepRibbon, RepRibbonNew, state);
1251 RepUpdateMacro(I, cRepCartoon, RepCartoonNew, state);
1252 RepUpdateMacro(I, cRepSurface, RepSurfaceNew, state);
1253 RepUpdateMacro(I, cRepLabel, RepLabelNew, state);
1254 RepUpdateMacro(I, cRepNonbonded, RepNonbondedNew, state);
1255 RepUpdateMacro(I, cRepNonbondedSphere, RepNonbondedSphereNew, state);
1256 RepUpdateMacro(I, cRepEllipsoid, RepEllipsoidNew, state);
1257
1258 for(a = 0; a < cRepCnt; a++)
1259 if(!I->Rep[a])
1260 I->Active[a] = false;
1261
1262 SceneInvalidate(G);
1263 OrthoBusyFast(G, 1, 1);
1264 if(Feedback(G, FB_CoordSet, FB_Blather)) {
1265 printf(" CoordSetUpdate-Leaving: object %s state %d cset %p\n",
1266 I->Obj->Name, state, (void *) I);
1267 }
1268 }
1269
1270
1271 /*========================================================================*/
CoordSetUpdateCoord2IdxMap(CoordSet * I,float cutoff)1272 void CoordSetUpdateCoord2IdxMap(CoordSet * I, float cutoff)
1273 {
1274 if(cutoff < R_SMALL4)
1275 cutoff = R_SMALL4;
1276 if(I->NIndex > 10) {
1277 if(I->Coord2Idx) {
1278 if((I->Coord2IdxDiv < cutoff) ||
1279 (((cutoff - I->Coord2IdxReq) / I->Coord2IdxReq) < -0.5F)) {
1280 MapFree(I->Coord2Idx);
1281 I->Coord2Idx = NULL;
1282 }
1283 }
1284 if(I->NIndex && (!I->Coord2Idx)) { /* NOTE: map based on stored coords */
1285 I->Coord2IdxReq = cutoff;
1286 I->Coord2IdxDiv = cutoff * 1.25F;
1287 I->Coord2Idx = MapNew(I->G, I->Coord2IdxDiv, I->Coord, I->NIndex, NULL);
1288 if(I->Coord2IdxDiv < I->Coord2Idx->Div)
1289 I->Coord2IdxDiv = I->Coord2Idx->Div;
1290 }
1291 }
1292 }
1293
1294 /*
1295 * RepToTransparencySetting: maps rep to transparency setting
1296 */
RepToTransparencySetting(const int rep_id)1297 static int RepToTransparencySetting(const int rep_id){
1298 switch(rep_id){
1299 // based on the rep,
1300 // transparency is set with different settings
1301 case cRepCyl:
1302 return cSetting_stick_transparency;
1303 case cRepSurface:
1304 return cSetting_transparency;
1305 case cRepSphere:
1306 return cSetting_sphere_transparency;
1307 case cRepEllipsoid:
1308 return cSetting_ellipsoid_transparency;
1309 case cRepCartoon:
1310 return cSetting_cartoon_transparency;
1311 case cRepRibbon:
1312 return cSetting_ribbon_transparency;
1313 }
1314 return 0;
1315 }
1316
1317 /*========================================================================*/
render(RenderInfo * info)1318 void CoordSet::render(RenderInfo * info)
1319 {
1320 CoordSet * I = this;
1321 PRINTFD(G, FB_CoordSet)
1322 " CoordSetRender: entered (%p).\n", (void *) I ENDFD;
1323
1324 if(!(info->ray || info->pick) &&
1325 (SettingGet_i(G, I->Setting, I->Obj->Setting,
1326 cSetting_defer_builds_mode) == 5)) {
1327 if(!info->pass) {
1328 ObjectUseColor((CObject *) I->Obj);
1329 if(I->Active[cRepLine])
1330 RepWireBondRenderImmediate(I, info);
1331 if(I->Active[cRepNonbonded])
1332 RepNonbondedRenderImmediate(I, info);
1333 if(I->Active[cRepSphere])
1334 RepSphereRenderImmediate(I, info);
1335 if(I->Active[cRepCyl])
1336 RepCylBondRenderImmediate(I, info);
1337 if(I->Active[cRepRibbon])
1338 RepRibbonRenderImmediate(I, info);
1339 }
1340 } else {
1341 int pass = info->pass;
1342 CRay *ray = info->ray;
1343 auto pick = info->pick;
1344 int a, aa, abit, aastart = 0, aaend = cRepCnt;
1345 ::Rep *r;
1346 int sculpt_vdw_vis_mode = SettingGet_i(G, I->Setting,
1347 I->Obj->Setting,
1348 cSetting_sculpt_vdw_vis_mode);
1349 if((!pass) && sculpt_vdw_vis_mode &&
1350 I->SculptCGO && (I->Obj->visRep & cRepCGOBit)) {
1351 if(ray) {
1352 int ok = CGORenderRay(I->SculptCGO, ray, info,
1353 ColorGet(G, I->Obj->Color), NULL, I->Setting, I->Obj->Setting);
1354 if (!ok){
1355 CGOFree(I->SculptCGO);
1356 CGOFree(I->SculptShaderCGO);
1357 }
1358 } else if(G->HaveGUI && G->ValidContext) {
1359 if(!pick) {
1360 int use_shader = SettingGetGlobal_b(G, cSetting_use_shaders);
1361 if (use_shader){
1362 if (!I->SculptShaderCGO){
1363 CGO *convertcgo = NULL;
1364 convertcgo = CGOCombineBeginEnd(I->SculptCGO, 0);
1365 if (convertcgo){
1366 I->SculptShaderCGO = CGOOptimizeToVBONotIndexed(convertcgo, 0);
1367 I->SculptShaderCGO->use_shader = true;
1368 CGOFree(convertcgo);
1369 }
1370 }
1371 } else {
1372 CGOFree(I->SculptShaderCGO);
1373 }
1374 if (I->SculptShaderCGO){
1375 CGORenderGL(I->SculptShaderCGO, NULL,
1376 I->Setting, I->Obj->Setting, info, NULL);
1377 } else {
1378 CGORenderGL(I->SculptCGO, NULL,
1379 I->Setting, I->Obj->Setting, info, NULL);
1380 }
1381 }
1382 }
1383 }
1384
1385 if (pick){
1386 int pick_labels = SettingGet_i(G, I->Setting,
1387 I->Obj->Setting,
1388 cSetting_pick_labels);
1389 if (pick_labels == 2){ // only pick labels
1390 aastart = cRepLabel;
1391 aaend = aastart + 1;
1392 }
1393 }
1394 for(aa = aastart; aa < aaend; aa++) {
1395 if(aa == cRepSurface) { /* reorder */
1396 a = cRepCell;
1397 } else if(aa == cRepCell) {
1398 a = cRepSurface;
1399 } else {
1400 a = aa;
1401 }
1402
1403 abit = (1 << a);
1404
1405 if(I->Active[a] && I->Rep[a]) {
1406 r = I->Rep[a];
1407 if(!ray) {
1408 ObjectUseColor((CObject *) I->Obj);
1409 } else {
1410 if(I->Obj)
1411 ray->wobble(
1412 SettingGet_i(G, I->Setting,
1413 I->Obj->Setting,
1414 cSetting_ray_texture),
1415 SettingGet_3fv(G, I->Setting,
1416 I->Obj->Setting,
1417 cSetting_ray_texture_settings));
1418 else
1419 ray->wobble(
1420 SettingGet_i(G, I->Setting,
1421 NULL, cSetting_ray_texture),
1422 SettingGet_3fv(G, I->Setting, NULL,
1423 cSetting_ray_texture_settings));
1424 ray->color3fv(ColorGet(G, I->Obj->Color));
1425 }
1426
1427 if(r->fRender) { /* do OpenGL rendering in three passes */
1428 if(ray || pick) {
1429
1430 /* here we need to iterate through and apply coordinate set matrices */
1431
1432 r->fRender(r, info);
1433 } else {
1434 bool t_mode_3 = SettingGetGlobal_i(G, cSetting_transparency_mode) == 3;
1435 bool render_both = t_mode_3; // render both opaque (1) and transparent (-1) pass if t_mode_3
1436 /* here we need to iterate through and apply coordinate set matrices */
1437
1438 switch (a) {
1439 case cRepVolume:
1440 {
1441 int t_mode = SettingGetGlobal_i(G, cSetting_transparency_mode);
1442 if (t_mode == 3 && pass == -1){
1443 r->fRender(r, info);
1444 }
1445 }
1446 break;
1447 case cRepLabel:
1448 {
1449 int t_mode_3 = SettingGetGlobal_i(G, cSetting_transparency_mode) == 3;
1450 if (pass == -1 || (t_mode_3 && pass == 1 /* TODO_OPENVR 0 */))
1451 r->fRender(r, info);
1452 }
1453 break;
1454 case cRepDot:
1455 case cRepCGO:
1456 case cRepCallback:
1457 if(pass == 1)
1458 r->fRender(r, info);
1459 break;
1460 case cRepLine:
1461 case cRepMesh:
1462 case cRepDash:
1463 case cRepCell:
1464 case cRepExtent:
1465 if(!pass)
1466 r->fRender(r, info);
1467 break;
1468 case cRepNonbonded:
1469 case cRepRibbon:
1470 render_both = false; // cartoon and nonbonded do not have atom-level transparency
1471 case cRepNonbondedSphere:
1472 case cRepSurface:
1473 case cRepEllipsoid:
1474 case cRepCyl:
1475 case cRepCartoon:
1476 case cRepSphere:
1477 {
1478 if (render_both){
1479 // for transparency_mode 3, render both opaque and transparent pass
1480 if (pass != 0){
1481 r->fRender(r, info);
1482 }
1483 } else {
1484 bool checkAlphaCGO = abit & (cRepSurfaceBit);
1485 int check_setting = RepToTransparencySetting(a);
1486 bool cont = true;
1487 if (checkAlphaCGO){
1488 if(info->alpha_cgo) {
1489 if(pass == 1){
1490 r->fRender(r, info);
1491 }
1492 cont = false;
1493 }
1494 }
1495 if (cont){
1496 if(check_setting && SettingGet_f(G, r->cs->Setting,
1497 r->obj->Setting, check_setting) > 0.0001) {
1498 /* if object has transparency, only render in -1 pass */
1499 if(pass == -1)
1500 r->fRender(r, info);
1501 } else if(pass == 1){
1502 r->fRender(r, info);
1503 }
1504 }
1505 }
1506 }
1507 break;
1508 }
1509 }
1510 }
1511 /* if(ray)
1512 ray->fWobble(ray,0,NULL); */
1513 }
1514 }
1515 }
1516 PRINTFD(G, FB_CoordSet)
1517 " CoordSetRender: leaving...\n" ENDFD;
1518 }
1519
1520
1521 /*========================================================================*/
CoordSet(PyMOLGlobals * G)1522 CoordSet::CoordSet(PyMOLGlobals* G)
1523 : CObjectState(G)
1524 {
1525 }
1526
1527 /*========================================================================*/
CoordSet(const CoordSet & cs)1528 CoordSet::CoordSet(const CoordSet& cs)
1529 : CObjectState(cs)
1530 {
1531 this->Obj = cs.Obj;
1532 this->Coord = cs.Coord;
1533 this->IdxToAtm = cs.IdxToAtm;
1534 this->NIndex = cs.NIndex;
1535 this->NAtIndex = cs.NAtIndex;
1536 this->prevNIndex = cs.prevNIndex;
1537 this->prevNAtIndex = cs.prevNAtIndex;
1538 std::copy(std::begin(cs.Rep), std::end(cs.Rep), std::begin(this->Rep));
1539 std::copy(std::begin(cs.Active), std::end(cs.Active), std::begin(this->Active));
1540 this->NTmpBond = cs.NTmpBond;
1541 this->NTmpLinkBond = cs.NTmpLinkBond;
1542 /* deep copy & return ptr to new symmetry */
1543 if(cs.Symmetry) {
1544 this->Symmetry = pymol::make_unique<CSymmetry>(*cs.Symmetry);
1545 }
1546 if(cs.PeriodicBox) {
1547 this->PeriodicBox = pymol::make_unique<CCrystal>(*cs.PeriodicBox);
1548 }
1549 std::copy(std::begin(cs.Name), std::end(cs.Name), std::begin(this->Name));
1550 this->PeriodicBoxType = cs.PeriodicBoxType;
1551 this->tmp_index = cs.tmp_index;
1552 this->Coord2IdxReq = cs.Coord2IdxReq;
1553 this->Coord2IdxDiv = cs.Coord2IdxDiv;
1554 this->objMolOpInvalidated = cs.objMolOpInvalidated;
1555
1556 // copy VLAs
1557 this->LabPos = cs.LabPos;
1558 this->RefPos = cs.RefPos;
1559 this->AtmToIdx = cs.AtmToIdx;
1560
1561 UtilZeroMem(this->Rep, sizeof(::Rep *) * cRepCnt);
1562
1563 #ifdef _PYMOL_IP_PROPERTIES
1564 #endif
1565
1566 }
1567
1568
1569 /*========================================================================*/
extendIndices(int nAtom)1570 int CoordSet::extendIndices(int nAtom)
1571 {
1572 CoordSet * I = this;
1573 int a, b;
1574 ObjectMolecule *obj = I->Obj;
1575 int ok = true;
1576 if(obj->DiscreteFlag) {
1577 ok = obj->setNDiscrete(nAtom);
1578
1579 if(I->AtmToIdx) { /* convert to discrete if necessary */
1580 I->AtmToIdx.freeP();
1581 if (ok){
1582 for(a = 0; a < I->NIndex; a++) {
1583 b = I->IdxToAtm[a];
1584 obj->DiscreteAtmToIdx[b] = a;
1585 obj->DiscreteCSet[b] = I;
1586 }
1587 }
1588 }
1589 }
1590 if(ok && I->NAtIndex < nAtom) {
1591 if(I->AtmToIdx) {
1592 VLASize(I->AtmToIdx, int, nAtom);
1593 CHECKOK(ok, I->AtmToIdx);
1594 if(ok && nAtom) {
1595 for(a = I->NAtIndex; a < nAtom; a++)
1596 I->AtmToIdx[a] = -1;
1597 }
1598 I->NAtIndex = nAtom;
1599 } else if(!obj->DiscreteFlag) {
1600 I->AtmToIdx = pymol::vla<int>(nAtom);
1601 CHECKOK(ok, I->AtmToIdx);
1602 if (ok){
1603 for(a = 0; a < nAtom; a++)
1604 I->AtmToIdx[a] = -1;
1605 I->NAtIndex = nAtom;
1606 }
1607 }
1608 }
1609 return ok;
1610 }
1611
1612
1613 /*========================================================================*/
appendIndices(int offset)1614 void CoordSet::appendIndices(int offset)
1615 {
1616 CoordSet * I = this;
1617 int a, b;
1618 ObjectMolecule *obj = I->Obj;
1619
1620 I->IdxToAtm = pymol::vla<int>(I->NIndex);
1621 if(I->NIndex) {
1622 ErrChkPtr(G, I->IdxToAtm);
1623 for(a = 0; a < I->NIndex; a++)
1624 I->IdxToAtm[a] = a + offset;
1625 }
1626 if(obj->DiscreteFlag) {
1627 VLACheck(obj->DiscreteAtmToIdx, int, I->NIndex + offset);
1628 VLACheck(obj->DiscreteCSet, CoordSet *, I->NIndex + offset);
1629 for(a = 0; a < I->NIndex; a++) {
1630 b = a + offset;
1631 obj->DiscreteAtmToIdx[b] = a;
1632 obj->DiscreteCSet[b] = I;
1633 }
1634 } else {
1635 I->AtmToIdx = pymol::vla<int>(I->NIndex + offset);
1636 if(I->NIndex + offset) {
1637 ErrChkPtr(G, I->AtmToIdx);
1638 for(a = 0; a < offset; a++)
1639 I->AtmToIdx[a] = -1;
1640 for(a = 0; a < I->NIndex; a++)
1641 I->AtmToIdx[a + offset] = a;
1642 }
1643 }
1644 I->NAtIndex = I->NIndex + offset;
1645 }
1646
1647
1648 /*========================================================================*/
enumIndices()1649 void CoordSet::enumIndices()
1650 {
1651 CoordSet * I = this;
1652 /* set up for simple case where 1 = 1, etc. */
1653 int a;
1654 I->AtmToIdx = pymol::vla<int>(I->NIndex);
1655 I->IdxToAtm = pymol::vla<int>(I->NIndex);
1656 if(I->NIndex) {
1657 ErrChkPtr(G, I->AtmToIdx);
1658 ErrChkPtr(G, I->IdxToAtm);
1659 for(a = 0; a < I->NIndex; a++) {
1660 I->AtmToIdx[a] = a;
1661 I->IdxToAtm[a] = a;
1662 }
1663 }
1664 I->NAtIndex = I->NIndex;
1665 }
1666
1667
1668 /*========================================================================*/
~CoordSet()1669 CoordSet::~CoordSet()
1670 {
1671 CoordSet * I = this;
1672 int a;
1673 ObjectMolecule *obj;
1674 if(I) {
1675 #ifdef _PYMOL_IP_PROPERTIES
1676 #endif
1677
1678 if (I->has_atom_state_settings){
1679 for(a = 0; a < I->NIndex; a++){
1680 if (I->has_atom_state_settings[a]){
1681 SettingUniqueDetachChain(G, I->atom_state_setting_id[a]);
1682 }
1683 }
1684 }
1685 for(a = 0; a < cRepCnt; a++)
1686 if(I->Rep[a])
1687 I->Rep[a]->fFree(I->Rep[a]);
1688 obj = I->Obj;
1689 if(obj)
1690 if(obj->DiscreteFlag) /* remove references to the atoms in discrete objects */
1691 for(a = 0; a < I->NIndex; a++) {
1692 obj->DiscreteAtmToIdx[I->IdxToAtm[a]] = -1;
1693 obj->DiscreteCSet[I->IdxToAtm[a]] = NULL;
1694 }
1695 MapFree(I->Coord2Idx);
1696 SettingFreeP(I->Setting);
1697 CGOFree(I->SculptCGO);
1698 }
1699 }
fFree()1700 void CoordSet::fFree()
1701 {
1702 delete this;
1703 }
LabPosTypeCopy(const LabPosType * src,LabPosType * dst)1704 void LabPosTypeCopy(const LabPosType * src, LabPosType * dst){
1705 dst->mode = src->mode;
1706 copy3f(src->pos, dst->pos);
1707 copy3f(src->offset, dst->offset);
1708 }
RefPosTypeCopy(const RefPosType * src,RefPosType * dst)1709 void RefPosTypeCopy(const RefPosType * src, RefPosType * dst){
1710 copy3f(src->coord, dst->coord);
1711 dst->specified = src->specified;
1712 }
1713
1714 #ifndef _PYMOL_NOPY
CoordSetSetSettingFromPyObject(PyMOLGlobals * G,CoordSet * cs,int at,int setting_id,PyObject * val)1715 int CoordSetSetSettingFromPyObject(PyMOLGlobals * G, CoordSet *cs, int at, int setting_id, PyObject *val){
1716 if (val == Py_None)
1717 val = NULL;
1718
1719 if (!val) {
1720 if (!cs->has_atom_state_settings ||
1721 !cs->has_atom_state_settings[at])
1722 return true;
1723 }
1724
1725 CoordSetCheckUniqueID(G, cs, at);
1726 cs->has_atom_state_settings[at] = true;
1727
1728 return SettingUniqueSetPyObject(G, cs->atom_state_setting_id[at], setting_id, val);
1729 }
1730 #endif
1731
CoordSetCheckSetting(PyMOLGlobals * G,CoordSet * cs,int at,int setting_id)1732 int CoordSetCheckSetting(PyMOLGlobals * G, CoordSet *cs, int at, int setting_id){
1733 if(!cs->has_atom_state_settings || !cs->has_atom_state_settings[at]) {
1734 return 0;
1735 } else {
1736 if(!SettingUniqueCheck(G, cs->atom_state_setting_id[at], setting_id)) {
1737 return 0;
1738 } else {
1739 return 1;
1740 }
1741 }
1742 }
1743
SettingGetIfDefinedPyObject(PyMOLGlobals * G,CoordSet * cs,int at,int setting_id)1744 PyObject *SettingGetIfDefinedPyObject(PyMOLGlobals * G, CoordSet *cs, int at, int setting_id){
1745 if(cs->has_atom_state_settings && cs->has_atom_state_settings[at]){
1746 return SettingUniqueGetPyObject(G, cs->atom_state_setting_id[at], setting_id);
1747 }
1748 return NULL;
1749 }
1750
CoordSetCheckUniqueID(PyMOLGlobals * G,CoordSet * I,int at)1751 int CoordSetCheckUniqueID(PyMOLGlobals * G, CoordSet *I, int at){
1752 if (!I->atom_state_setting_id){
1753 I->atom_state_setting_id = pymol::vla<int>(I->NIndex);
1754 }
1755 if (!I->has_atom_state_settings){
1756 I->has_atom_state_settings = pymol::vla<char>(I->NIndex);
1757 }
1758 if (!I->atom_state_setting_id[at]){
1759 I->atom_state_setting_id[at] = AtomInfoGetNewUniqueID(G);
1760 }
1761 return I->atom_state_setting_id[at];
1762 }
1763
1764 template <typename V>
AtomStateGetSetting(ATOMSTATEGETSETTINGARGS,V * out)1765 void AtomStateGetSetting(ATOMSTATEGETSETTINGARGS, V * out) {
1766 if (cs->has_atom_state_settings &&
1767 cs->has_atom_state_settings[idx] &&
1768 SettingUniqueGetIfDefined(G, cs->atom_state_setting_id[idx], setting_id, out))
1769 return;
1770
1771 if (AtomSettingGetIfDefined(G, ai, setting_id, out))
1772 return;
1773
1774 *out = SettingGet<V>(G, cs->Setting, obj->Setting, setting_id);
1775 }
1776
1777 template void AtomStateGetSetting(ATOMSTATEGETSETTINGARGS, int * out);
1778 template void AtomStateGetSetting(ATOMSTATEGETSETTINGARGS, float * out);
1779 template void AtomStateGetSetting(ATOMSTATEGETSETTINGARGS, const float ** out);
1780