Context Navigation

source: vtkvis/trunk/Molecule.cpp @ 6226

Last change on this file since 6226 was 5835, checked in by ldelgass, 9 years ago
Require VTK >= 6.0.0
Property svn:eol-style set to `native`
File size: 34.2 KB

Line
1	/* -- mode: c++; c-basic-offset: 4; indent-tabs-mode: nil -- */
2	/*
3	* Copyright (C) 2004-2012 HUBzero Foundation, LLC
4	*
5	* Author: Leif Delgass <ldelgass@purdue.edu>
6	*/
7
8	#include <cstdio>
9	#include <cfloat>
10	#include <cassert>
11
12	#include <vtkDataSet.h>
13	#include <vtkCellArray.h>
14	#include <vtkPointData.h>
15	#include <vtkFloatArray.h>
16	#include <vtkDoubleArray.h>
17	#include <vtkIntArray.h>
18	#include <vtkStringArray.h>
19	#include <vtkPolyData.h>
20	#include <vtkPolyDataMapper.h>
21	#include <vtkActor.h>
22	#include <vtkProperty.h>
23	#include <vtkSphereSource.h>
24	#include <vtkCylinderSource.h>
25	#include <vtkTransform.h>
26	#include <vtkTransformPolyDataFilter.h>
27	#include <vtkGlyph3DMapper.h>
28	#include <vtkPointSetToLabelHierarchy.h>
29	#include <vtkLabelPlacementMapper.h>
30	#include <vtkTextProperty.h>
31	#include <vtkTransformPolyDataFilter.h>
32
33	#include "Molecule.h"
34	#include "MoleculeData.h"
35	#include "Renderer.h"
36	#include "Trace.h"
37
38	using namespace VtkVis;
39
40	Molecule::Molecule() :
41	GraphicsObject(),
42	_radiusScale(0.3),
43	_atomScaling(COVALENT_RADIUS),
44	_labelsOn(false),
45	_colorMap(NULL),
46	_colorMode(COLOR_BY_ELEMENTS),
47	_colorFieldType(DataSet::POINT_DATA)
48	{
49	_bondColor[0] = _bondColor[1] = _bondColor[2] = 1.0f;
50	_colorFieldRange[0] = DBL_MAX;
51	_colorFieldRange[1] = -DBL_MAX;
52	_faceCulling = true;
53	}
54
55	Molecule::~Molecule()
56	{
57	#ifdef WANT_TRACE
58	if (_dataSet != NULL)
59	TRACE("Deleting Molecule for %s", _dataSet->getName().c_str());
60	else
61	TRACE("Deleting Molecule with NULL DataSet");
62	#endif
63	}
64
65	void Molecule::setDataSet(DataSet *dataSet,
66	Renderer *renderer)
67	{
68	if (_dataSet != dataSet) {
69	_dataSet = dataSet;
70
71	_renderer = renderer;
72
73	if (renderer->getUseCumulativeRange()) {
74	renderer->getCumulativeDataRange(_dataRange,
75	_dataSet->getActiveScalarsName(),
76	1);
77	renderer->getCumulativeDataRange(_vectorMagnitudeRange,
78	_dataSet->getActiveVectorsName(),
79	3);
80	for (int i = 0; i < 3; i++) {
81	renderer->getCumulativeDataRange(_vectorComponentRange[i],
82	_dataSet->getActiveVectorsName(),
83	3, i);
84	}
85	} else {
86	_dataSet->getScalarRange(_dataRange);
87	_dataSet->getVectorRange(_vectorMagnitudeRange);
88	for (int i = 0; i < 3; i++) {
89	_dataSet->getVectorRange(_vectorComponentRange[i], i);
90	}
91	}
92
93	update();
94	}
95	}
96
97	/**
98	* \brief Create and initialize VTK Props to render a Molecule
99	*/
100	void Molecule::initProp()
101	{
102	if (_atomProp == NULL) {
103	_atomProp = vtkSmartPointer<vtkActor>::New();
104	if (_faceCulling && _opacity == 1.0)
105	setCulling(_atomProp->GetProperty(), true);
106	_atomProp->GetProperty()->EdgeVisibilityOff();
107	_atomProp->GetProperty()->SetEdgeColor(_edgeColor[0], _edgeColor[1], _edgeColor[2]);
108	_atomProp->GetProperty()->SetLineWidth(_edgeWidth);
109	_atomProp->GetProperty()->SetOpacity(_opacity);
110	_atomProp->GetProperty()->SetAmbient(.2);
111	_atomProp->GetProperty()->SetSpecular(.2);
112	_atomProp->GetProperty()->SetSpecularPower(80.0);
113	if (!_lighting)
114	_atomProp->GetProperty()->LightingOff();
115	}
116	if (_bondProp == NULL) {
117	_bondProp = vtkSmartPointer<vtkActor>::New();
118	if (_faceCulling && _opacity == 1.0)
119	setCulling(_bondProp->GetProperty(), true);
120	_bondProp->GetProperty()->EdgeVisibilityOff();
121	_bondProp->GetProperty()->SetColor(_bondColor[0], _bondColor[1], _bondColor[2]);
122	_bondProp->GetProperty()->SetEdgeColor(_edgeColor[0], _edgeColor[1], _edgeColor[2]);
123	_bondProp->GetProperty()->SetLineWidth(_edgeWidth);
124	_bondProp->GetProperty()->SetOpacity(_opacity);
125	_bondProp->GetProperty()->SetAmbient(.2);
126	_bondProp->GetProperty()->SetSpecular(.2);
127	_bondProp->GetProperty()->SetSpecularPower(80.0);
128	if (!_lighting)
129	_bondProp->GetProperty()->LightingOff();
130	}
131	if (_labelProp == NULL) {
132	_labelProp = vtkSmartPointer<vtkActor2D>::New();
133	}
134	if (_prop == NULL) {
135	_prop = vtkSmartPointer<vtkAssembly>::New();
136	}
137	}
138
139	/**
140	* \brief Internal method to set up pipeline after a state change
141	*/
142	void Molecule::update()
143	{
144	if (_dataSet == NULL) {
145	return;
146	}
147	vtkDataSet *ds = _dataSet->getVtkDataSet();
148
149	if (_atomMapper == NULL) {
150	_atomMapper = vtkSmartPointer<vtkGlyph3DMapper>::New();
151	_atomMapper->SetResolveCoincidentTopologyToPolygonOffset();
152	_atomMapper->ScalarVisibilityOn();
153	}
154	if (_bondMapper == NULL) {
155	_bondMapper = vtkSmartPointer<vtkGlyph3DMapper>::New();
156	_bondMapper->SetResolveCoincidentTopologyToPolygonOffset();
157	_bondMapper->ScalarVisibilityOn();
158	}
159	if (_labelMapper == NULL) {
160	_labelMapper = vtkSmartPointer<vtkLabelPlacementMapper>::New();
161	_labelMapper->SetShapeToRoundedRect();
162	_labelMapper->SetBackgroundColor(1.0, 1.0, 0.7);
163	_labelMapper->SetBackgroundOpacity(0.8);
164	_labelMapper->SetMargin(3);
165	}
166
167	if (_lut == NULL) {
168	if (ds->GetPointData() == NULL \|\|
169	ds->GetPointData()->GetScalars() == NULL \|\|
170	strcmp(ds->GetPointData()->GetScalars()->GetName(), "element") != 0) {
171	WARN("No element array in dataset %s", _dataSet->getName().c_str());
172	setColorMap(ColorMap::getDefault());
173	if (_atomScaling != NO_ATOM_SCALING)
174	_atomScaling = NO_ATOM_SCALING;
175	} else {
176	TRACE("Using element default colormap");
177	setColorMap(ColorMap::getElementDefault());
178	}
179	}
180
181	initProp();
182
183	addLabelArray(ds);
184	addRadiusArray(ds, _atomScaling, _radiusScale);
185	computeBonds(ds);
186
187	vtkPolyData *pd = vtkPolyData::SafeDownCast(ds);
188	if (pd) {
189	TRACE("Points: %d Verts: %d Lines: %d Polys: %d Strips: %d",
190	pd->GetNumberOfPoints(),
191	pd->GetNumberOfVerts(),
192	pd->GetNumberOfLines(),
193	pd->GetNumberOfPolys(),
194	pd->GetNumberOfStrips());
195	// DataSet is a vtkPolyData
196	if (pd->GetNumberOfLines() > 0) {
197	// Bonds
198	setupBondPolyData();
199
200	if (_cylinderSource == NULL) {
201	_cylinderSource = vtkSmartPointer<vtkCylinderSource>::New();
202	_cylinderSource->SetRadius(0.075);
203	_cylinderSource->SetHeight(1.0);
204	_cylinderSource->SetResolution(12);
205	_cylinderSource->CappingOff();
206	_cylinderTrans = vtkSmartPointer<vtkTransformPolyDataFilter>::New();
207	_cylinderTrans->SetInputConnection(_cylinderSource->GetOutputPort());
208	vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
209	trans->RotateZ(-90.0);
210	_cylinderTrans->SetTransform(trans);
211	}
212	if (_lineSource == NULL) {
213	_lineSource = vtkSmartPointer<vtkLineSource>::New();
214	_lineSource->SetPoint1(-0.5, 0, 0);
215	_lineSource->SetPoint2(0.5, 0, 0);
216	}
217
218	_bondMapper->SetSourceConnection(_cylinderTrans->GetOutputPort());
219	_bondMapper->SetInputData(_bondPD);
220	_bondMapper->SetOrientationArray("bond_orientations");
221	_bondMapper->SetOrientationModeToDirection();
222	_bondMapper->OrientOn();
223	_bondMapper->SetScaleArray("bond_scales");
224	_bondMapper->SetScaleModeToScaleByVectorComponents();
225	_bondMapper->ScalingOn();
226	_bondMapper->ClampingOff();
227
228	_bondProp->SetMapper(_bondMapper);
229	getAssembly()->AddPart(_bondProp);
230	}
231	if (pd->GetNumberOfPoints() > 0) {
232	if (_labelHierarchy == NULL) {
233	_labelHierarchy = vtkSmartPointer<vtkPointSetToLabelHierarchy>::New();
234	}
235	if (_labelTransform == NULL) {
236	_labelTransform = vtkSmartPointer<vtkTransform>::New();
237	}
238	vtkSmartPointer<vtkTransformPolyDataFilter> transformFilter = vtkSmartPointer<vtkTransformPolyDataFilter>::New();
239	transformFilter->SetInputData(pd);
240	transformFilter->SetTransform(_labelTransform);
241	_labelHierarchy->SetInputConnection(transformFilter->GetOutputPort());
242	_labelHierarchy->SetLabelArrayName("_atom_labels");
243	_labelHierarchy->GetTextProperty()->SetColor(0, 0, 0);
244	_labelMapper->SetInputConnection(_labelHierarchy->GetOutputPort());
245	_labelProp->SetMapper(_labelMapper);
246
247	// Atoms
248	if (_sphereSource == NULL) {
249	_sphereSource = vtkSmartPointer<vtkSphereSource>::New();
250	_sphereSource->SetRadius(1.0);
251	_sphereSource->SetThetaResolution(14);
252	_sphereSource->SetPhiResolution(14);
253	}
254
255	_atomMapper->SetSourceConnection(_sphereSource->GetOutputPort());
256	_atomMapper->SetInputData(pd);
257	_atomMapper->SetScaleArray("_radii");
258	_atomMapper->SetScaleModeToScaleByMagnitude();
259	_atomMapper->ScalingOn();
260	_atomMapper->OrientOff();
261
262	_atomProp->SetMapper(_atomMapper);
263	getAssembly()->AddPart(_atomProp);
264	}
265	} else {
266	// DataSet is NOT a vtkPolyData
267	ERROR("DataSet is not a PolyData");
268	return;
269	}
270
271	setAtomLabelVisibility(_labelsOn);
272
273	if (pd->GetNumberOfPoints() > 0) {
274	_atomMapper->Update();
275	_labelMapper->Update();
276	}
277	if (pd->GetNumberOfLines() > 0) {
278	_bondMapper->Update();
279	}
280	}
281
282	void Molecule::updateLabelTransform()
283	{
284	if (_labelTransform != NULL && getAssembly() != NULL) {
285	_labelTransform->SetMatrix(getAssembly()->GetMatrix());
286	}
287	}
288
289	void Molecule::updateRanges(Renderer *renderer)
290	{
291	if (_dataSet == NULL) {
292	ERROR("called before setDataSet");
293	return;
294	}
295
296	if (renderer->getUseCumulativeRange()) {
297	renderer->getCumulativeDataRange(_dataRange,
298	_dataSet->getActiveScalarsName(),
299	1);
300	renderer->getCumulativeDataRange(_vectorMagnitudeRange,
301	_dataSet->getActiveVectorsName(),
302	3);
303	for (int i = 0; i < 3; i++) {
304	renderer->getCumulativeDataRange(_vectorComponentRange[i],
305	_dataSet->getActiveVectorsName(),
306	3, i);
307	}
308	} else {
309	_dataSet->getScalarRange(_dataRange);
310	_dataSet->getVectorRange(_vectorMagnitudeRange);
311	for (int i = 0; i < 3; i++) {
312	_dataSet->getVectorRange(_vectorComponentRange[i], i);
313	}
314	}
315
316	// Need to update color map ranges
317	double *rangePtr = _colorFieldRange;
318	if (_colorFieldRange[0] > _colorFieldRange[1]) {
319	rangePtr = NULL;
320	}
321	setColorMode(_colorMode, _colorFieldType, _colorFieldName.c_str(), rangePtr);
322	}
323
324	void Molecule::setColorMode(ColorMode mode)
325	{
326	_colorMode = mode;
327	if (_dataSet == NULL)
328	return;
329
330	switch (mode) {
331	case COLOR_BY_ELEMENTS: // Assume default scalar is "element" array
332	setColorMode(mode,
333	DataSet::POINT_DATA,
334	"element");
335	break;
336	case COLOR_BY_SCALAR:
337	setColorMode(mode,
338	_dataSet->getActiveScalarsType(),
339	_dataSet->getActiveScalarsName());
340	break;
341	case COLOR_BY_VECTOR_MAGNITUDE:
342	setColorMode(mode,
343	_dataSet->getActiveVectorsType(),
344	_dataSet->getActiveVectorsName());
345	break;
346	case COLOR_BY_VECTOR_X:
347	setColorMode(mode,
348	_dataSet->getActiveVectorsType(),
349	_dataSet->getActiveVectorsName());
350	break;
351	case COLOR_BY_VECTOR_Y:
352	setColorMode(mode,
353	_dataSet->getActiveVectorsType(),
354	_dataSet->getActiveVectorsName());
355	break;
356	case COLOR_BY_VECTOR_Z:
357	setColorMode(mode,
358	_dataSet->getActiveVectorsType(),
359	_dataSet->getActiveVectorsName());
360	break;
361	case COLOR_CONSTANT:
362	default:
363	setColorMode(mode, DataSet::POINT_DATA, NULL, NULL);
364	break;
365	}
366	}
367
368	void Molecule::setColorMode(ColorMode mode,
369	const char *name, double range[2])
370	{
371	if (_dataSet == NULL)
372	return;
373	DataSet::DataAttributeType type = DataSet::POINT_DATA;
374	int numComponents = 1;
375	if (name != NULL && strlen(name) > 0 &&
376	!_dataSet->getFieldInfo(name, &type, &numComponents)) {
377	ERROR("Field not found: %s", name);
378	return;
379	}
380	setColorMode(mode, type, name, range);
381	}
382
383	void Molecule::setColorMode(ColorMode mode, DataSet::DataAttributeType type,
384	const char *name, double range[2])
385	{
386	_colorMode = mode;
387	_colorFieldType = type;
388	if (name == NULL)
389	_colorFieldName.clear();
390	else
391	_colorFieldName = name;
392	if (range == NULL) {
393	_colorFieldRange[0] = DBL_MAX;
394	_colorFieldRange[1] = -DBL_MAX;
395	} else {
396	memcpy(_colorFieldRange, range, sizeof(double)*2);
397	}
398
399	if (_dataSet == NULL \|\| (_atomMapper == NULL && _bondMapper == NULL))
400	return;
401
402	switch (type) {
403	case DataSet::POINT_DATA:
404	if (_atomMapper != NULL)
405	_atomMapper->SetScalarModeToUsePointFieldData();
406	if (_bondMapper != NULL)
407	_bondMapper->SetScalarModeToUsePointFieldData();
408	break;
409	case DataSet::CELL_DATA:
410	if (_atomMapper != NULL)
411	_atomMapper->SetScalarModeToUseCellFieldData();
412	if (_bondMapper != NULL)
413	_bondMapper->SetScalarModeToUseCellFieldData();
414	break;
415	default:
416	ERROR("Unsupported DataAttributeType: %d", type);
417	return;
418	}
419
420	if (name != NULL && strlen(name) > 0) {
421	if (_atomMapper != NULL)
422	_atomMapper->SelectColorArray(name);
423	if (_bondMapper != NULL)
424	_bondMapper->SelectColorArray(name);
425	} else {
426	if (_atomMapper != NULL)
427	_atomMapper->SetScalarModeToDefault();
428	if (_bondMapper != NULL)
429	_bondMapper->SetScalarModeToDefault();
430	}
431
432	if (_lut != NULL) {
433	if (range != NULL) {
434	_lut->SetRange(range);
435	} else if (mode == COLOR_BY_ELEMENTS) {
436	double range[2];
437	range[0] = 0;
438	range[1] = NUM_ELEMENTS;
439	_lut->SetRange(range);
440	} else if (name != NULL && strlen(name) > 0) {
441	double r[2];
442	int comp = -1;
443	if (mode == COLOR_BY_VECTOR_X)
444	comp = 0;
445	else if (mode == COLOR_BY_VECTOR_Y)
446	comp = 1;
447	else if (mode == COLOR_BY_VECTOR_Z)
448	comp = 2;
449
450	if (_renderer->getUseCumulativeRange()) {
451	int numComponents;
452	if (!_dataSet->getFieldInfo(name, type, &numComponents)) {
453	ERROR("Field not found: %s, type: %d", name, type);
454	return;
455	} else if (numComponents < comp+1) {
456	ERROR("Request for component %d in field with %d components",
457	comp, numComponents);
458	return;
459	}
460	_renderer->getCumulativeDataRange(r, name, type, numComponents, comp);
461	} else {
462	_dataSet->getDataRange(r, name, type, comp);
463	}
464	_lut->SetRange(r);
465	}
466	}
467
468	switch (mode) {
469	case COLOR_BY_ELEMENTS:
470	if (_atomMapper != NULL)
471	_atomMapper->ScalarVisibilityOn();
472	if (_bondMapper != NULL)
473	_bondMapper->ScalarVisibilityOn();
474	break;
475	case COLOR_BY_SCALAR:
476	if (_atomMapper != NULL)
477	_atomMapper->ScalarVisibilityOn();
478	if (_bondMapper != NULL)
479	_bondMapper->ScalarVisibilityOn();
480	break;
481	case COLOR_BY_VECTOR_MAGNITUDE:
482	if (_atomMapper != NULL)
483	_atomMapper->ScalarVisibilityOn();
484	if (_bondMapper != NULL)
485	_bondMapper->ScalarVisibilityOn();
486	if (_lut != NULL) {
487	_lut->SetVectorModeToMagnitude();
488	}
489	break;
490	case COLOR_BY_VECTOR_X:
491	if (_atomMapper != NULL)
492	_atomMapper->ScalarVisibilityOn();
493	if (_bondMapper != NULL)
494	_bondMapper->ScalarVisibilityOn();
495	if (_lut != NULL) {
496	_lut->SetVectorModeToComponent();
497	_lut->SetVectorComponent(0);
498	}
499	break;
500	case COLOR_BY_VECTOR_Y:
501	if (_atomMapper != NULL)
502	_atomMapper->ScalarVisibilityOn();
503	if (_bondMapper != NULL)
504	_bondMapper->ScalarVisibilityOn();
505	if (_lut != NULL) {
506	_lut->SetVectorModeToComponent();
507	_lut->SetVectorComponent(1);
508	}
509	break;
510	case COLOR_BY_VECTOR_Z:
511	if (_atomMapper != NULL)
512	_atomMapper->ScalarVisibilityOn();
513	if (_bondMapper != NULL)
514	_bondMapper->ScalarVisibilityOn();
515	if (_lut != NULL) {
516	_lut->SetVectorModeToComponent();
517	_lut->SetVectorComponent(2);
518	}
519	break;
520	case COLOR_CONSTANT:
521	default:
522	if (_atomMapper != NULL)
523	_atomMapper->ScalarVisibilityOff();
524	if (_bondMapper != NULL)
525	_bondMapper->ScalarVisibilityOff();
526	break;
527	}
528	}
529
530	/**
531	* \brief Called when the color map has been edited
532	*/
533	void Molecule::updateColorMap()
534	{
535	setColorMap(_colorMap);
536	}
537
538	/**
539	* \brief Associate a colormap lookup table with the DataSet
540	*/
541	void Molecule::setColorMap(ColorMap *cmap)
542	{
543	if (cmap == NULL)
544	return;
545
546	_colorMap = cmap;
547
548	if (_lut == NULL) {
549	_lut = vtkSmartPointer<vtkLookupTable>::New();
550	if (_atomMapper != NULL) {
551	_atomMapper->UseLookupTableScalarRangeOn();
552	_atomMapper->SetLookupTable(_lut);
553	}
554	if (_bondMapper != NULL) {
555	_bondMapper->UseLookupTableScalarRangeOn();
556	_bondMapper->SetLookupTable(_lut);
557	}
558	_lut->DeepCopy(cmap->getLookupTable());
559	switch (_colorMode) {
560	case COLOR_BY_ELEMENTS: {
561	double range[2];
562	range[0] = 0;
563	range[1] = NUM_ELEMENTS;
564	_lut->SetRange(range);
565	}
566	break;
567	case COLOR_CONSTANT:
568	case COLOR_BY_SCALAR:
569	_lut->SetRange(_dataRange);
570	break;
571	case COLOR_BY_VECTOR_MAGNITUDE:
572	_lut->SetRange(_vectorMagnitudeRange);
573	break;
574	case COLOR_BY_VECTOR_X:
575	_lut->SetRange(_vectorComponentRange[0]);
576	break;
577	case COLOR_BY_VECTOR_Y:
578	_lut->SetRange(_vectorComponentRange[1]);
579	break;
580	case COLOR_BY_VECTOR_Z:
581	_lut->SetRange(_vectorComponentRange[2]);
582	break;
583	default:
584	break;
585	}
586	} else {
587	double range[2];
588	_lut->GetTableRange(range);
589	_lut->DeepCopy(cmap->getLookupTable());
590	_lut->SetRange(range);
591	_lut->Modified();
592	}
593	}
594
595	void Molecule::setAtomLabelField(const char *fieldName)
596	{
597	if (_labelHierarchy != NULL) {
598	if (strcmp(fieldName, "default") == 0) {
599	_labelHierarchy->SetLabelArrayName("_atom_labels");
600	} else {
601	_labelHierarchy->SetLabelArrayName(fieldName);
602	}
603	}
604	}
605
606	/**
607	* \brief Turn on/off rendering of atom labels
608	*/
609	void Molecule::setAtomLabelVisibility(bool state)
610	{
611	_labelsOn = state;
612	if (_labelProp != NULL) {
613	_labelProp->SetVisibility((state ? 1 : 0));
614	}
615	}
616
617	/**
618	* \brief Turn on/off rendering of the atoms
619	*/
620	void Molecule::setAtomVisibility(bool state)
621	{
622	if (_atomProp != NULL) {
623	_atomProp->SetVisibility((state ? 1 : 0));
624	}
625	}
626
627	/**
628	* \brief Turn on/off rendering of the bonds
629	*/
630	void Molecule::setBondVisibility(bool state)
631	{
632	if (_bondProp != NULL) {
633	_bondProp->SetVisibility((state ? 1 : 0));
634	}
635	}
636
637	/**
638	* \brief Toggle visibility of the prop
639	*/
640	void Molecule::setVisibility(bool state)
641	{
642	GraphicsObject::setVisibility(state);
643	if (_labelProp != NULL) {
644	if (!state)
645	_labelProp->SetVisibility(0);
646	else
647	setAtomLabelVisibility(_labelsOn);
648	}
649	}
650
651	/**
652	* \brief Set opacity of molecule
653	*/
654	void Molecule::setOpacity(double opacity)
655	{
656	GraphicsObject::setOpacity(opacity);
657	if (_labelMapper != NULL) {
658	_labelMapper->SetBackgroundOpacity(opacity);
659	}
660	}
661
662	void Molecule::setAtomQuality(double quality)
663	{
664	if (_sphereSource == NULL)
665	return;
666
667	if (quality > 10.0)
668	quality = 10.0;
669
670	int thetaRes = (int)(quality * 14.0);
671	int phiRes = (int)(quality * 14.0);
672	if (thetaRes < 4) thetaRes = 4;
673	if (phiRes < 3) phiRes = 3;
674
675	_sphereSource->SetThetaResolution(thetaRes);
676	_sphereSource->SetPhiResolution(phiRes);
677
678	if (_atomMapper != NULL) {
679	_atomMapper->Modified();
680	_atomMapper->Update();
681	}
682	}
683
684	void Molecule::setBondQuality(double quality)
685	{
686	if (_cylinderSource == NULL)
687	return;
688
689	if (quality > 10.0)
690	quality = 10.0;
691
692	int res = (int)(quality * 12.0);
693	if (res < 3) res = 3;
694
695	_cylinderSource->SetResolution(res);
696
697	if (_bondMapper != NULL) {
698	_bondMapper->Modified();
699	_bondMapper->Update();
700	}
701	}
702
703	void Molecule::setBondStyle(BondStyle style)
704	{
705	switch (style) {
706	case BOND_STYLE_CYLINDER:
707	if (_bondProp != NULL) {
708	_bondProp->GetProperty()->SetLineWidth(_edgeWidth);
709	_bondProp->GetProperty()->SetLighting(_lighting ? 1 : 0);
710	}
711	if (_bondMapper != NULL && _cylinderTrans != NULL &&
712	_bondMapper->GetInputConnection(1, 0) != _cylinderTrans->GetOutputPort()) {
713	_cylinderTrans->Modified();
714	_bondMapper->SetSourceConnection(_cylinderTrans->GetOutputPort());
715	_bondMapper->Modified();
716	}
717	break;
718	case BOND_STYLE_LINE:
719	if (_bondProp != NULL) {
720	_bondProp->GetProperty()->LightingOff();
721	}
722	if (_bondMapper != NULL && _lineSource != NULL &&
723	_bondMapper->GetInputConnection(1, 0) != _lineSource->GetOutputPort()) {
724	_lineSource->Modified();
725	_bondMapper->SetSourceConnection(_lineSource->GetOutputPort());
726	_bondMapper->Modified();
727	}
728	break;
729	default:
730	WARN("Unknown bond style");
731	}
732	}
733
734	/**
735	* \brief Set constant bond color
736	*/
737	void Molecule::setBondColor(float color[3])
738	{
739	_bondColor[0] = color[0];
740	_bondColor[1] = color[1];
741	_bondColor[2] = color[2];
742	if (_bondProp != NULL) {
743	_bondProp->GetProperty()->SetColor(_bondColor[0], _bondColor[1], _bondColor[2]);
744	}
745	}
746
747	/**
748	* \brief Set mode to determine how bonds are colored
749	*/
750	void Molecule::setBondColorMode(BondColorMode mode)
751	{
752	if (_bondMapper == NULL) return;
753
754	switch (mode) {
755	case BOND_COLOR_BY_ELEMENTS:
756	_bondMapper->ScalarVisibilityOn();
757	break;
758	case BOND_COLOR_CONSTANT:
759	_bondMapper->ScalarVisibilityOff();
760	break;
761	default:
762	WARN("Unknown bond color mode");
763	}
764	}
765
766	/**
767	* \brief Set a group of world coordinate planes to clip rendering
768	*
769	* Passing NULL for planes will remove all cliping planes
770	*/
771	void Molecule::setClippingPlanes(vtkPlaneCollection *planes)
772	{
773	if (_atomMapper != NULL) {
774	_atomMapper->SetClippingPlanes(planes);
775	}
776	if (_bondMapper != NULL) {
777	_bondMapper->SetClippingPlanes(planes);
778	}
779	}
780
781	/**
782	* \brief Set the radius type used for scaling atoms
783	*/
784	void Molecule::setAtomScaling(AtomScaling state)
785	{
786	if (state == _atomScaling)
787	return;
788
789	_atomScaling = state;
790	if (_dataSet != NULL) {
791	vtkDataSet *ds = _dataSet->getVtkDataSet();
792	addRadiusArray(ds, _atomScaling, _radiusScale);
793	if (_atomMapper != NULL) {
794	assert(ds->GetPointData() != NULL &&
795	ds->GetPointData()->GetArray("_radii") != NULL);
796	_atomMapper->SetScaleModeToScaleByMagnitude();
797	_atomMapper->SetScaleArray("_radii");
798	_atomMapper->ScalingOn();
799	}
800	}
801	}
802
803	/**
804	* \brief Set the constant radius scaling factor for atoms. This
805	* can be used to convert from Angstroms to atom coordinates units.
806	*/
807	void Molecule::setAtomRadiusScale(double scale)
808	{
809	if (scale == _radiusScale)
810	return;
811
812	_radiusScale = scale;
813	if (_dataSet != NULL) {
814	vtkDataSet *ds = _dataSet->getVtkDataSet();
815	addRadiusArray(ds, _atomScaling, _radiusScale);
816	}
817	}
818
819	/**
820	* \brief Set the constant radius scaling factor for bonds.
821	*/
822	void Molecule::setBondRadiusScale(double scale)
823	{
824	if (_cylinderSource != NULL &&
825	_cylinderSource->GetRadius() != scale) {
826	_cylinderSource->SetRadius(scale);
827	// Workaround bug with source modification not causing
828	// mapper to be updated
829	if (_bondMapper != NULL) {
830	_bondMapper->Modified();
831	}
832	}
833	}
834
835	void Molecule::setupBondPolyData()
836	{
837	if (_dataSet == NULL)
838	return;
839	if (_bondPD == NULL) {
840	_bondPD = vtkSmartPointer<vtkPolyData>::New();
841	} else {
842	_bondPD->Initialize();
843	}
844	vtkPolyData *pd = vtkPolyData::SafeDownCast(_dataSet->getVtkDataSet());
845	if (pd == NULL)
846	return;
847	vtkCellArray *lines = pd->GetLines();
848	lines->InitTraversal();
849	vtkIdType npts, *pts;
850	vtkSmartPointer<vtkPoints> bondPoints = vtkSmartPointer<vtkPoints>::New();
851	vtkSmartPointer<vtkIntArray> bondElements = vtkSmartPointer<vtkIntArray>::New();
852	vtkSmartPointer<vtkDoubleArray> bondVectors = vtkSmartPointer<vtkDoubleArray>::New();
853	vtkSmartPointer<vtkDoubleArray> bondScales = vtkSmartPointer<vtkDoubleArray>::New();
854	bondElements->SetName("element");
855	bondVectors->SetName("bond_orientations");
856	bondVectors->SetNumberOfComponents(3);
857	bondScales->SetName("bond_scales");
858	bondScales->SetNumberOfComponents(3);
859	vtkDataArray *elements = NULL;
860	if (pd->GetPointData() != NULL &&
861	pd->GetPointData()->GetScalars() != NULL &&
862	strcmp(pd->GetPointData()->GetScalars()->GetName(), "element") == 0) {
863	elements = pd->GetPointData()->GetScalars();
864	}
865	for (int i = 0; lines->GetNextCell(npts, pts); i++) {
866	assert(npts == 2);
867	double pt0[3], pt1[3];
868	double newPt0[3], newPt1[3];
869	//double *pt0 = pd->GetPoint(pts[0]);
870	//double *pt1 = pd->GetPoint(pts[1]);
871	pd->GetPoint(pts[0], pt0);
872	pd->GetPoint(pts[1], pt1);
873	double center[3];
874
875	for (int j = 0; j < 3; j++)
876	center[j] = pt0[j] + (pt1[j] - pt0[j]) * 0.5;
877	for (int j = 0; j < 3; j++)
878	newPt0[j] = pt0[j] + (pt1[j] - pt0[j]) * 0.25;
879	for (int j = 0; j < 3; j++)
880	newPt1[j] = pt0[j] + (pt1[j] - pt0[j]) * 0.75;
881
882	bondPoints->InsertNextPoint(newPt0);
883	bondPoints->InsertNextPoint(newPt1);
884
885	#ifdef DEBUG
886	TRACE("Bond %d: (%g,%g,%g)-(%g,%g,%g)-(%g,%g,%g)", i,
887	pt0[0], pt0[1], pt0[2],
888	center[0], center[1], center[2],
889	pt1[0], pt1[1], pt1[2]);
890	#endif
891
892	double vec[3];
893	for (int j = 0; j < 3; j++)
894	vec[j] = center[j] - pt0[j];
895
896	bondVectors->InsertNextTupleValue(vec);
897	bondVectors->InsertNextTupleValue(vec);
898	#ifdef DEBUG
899	TRACE("Bond %d, vec: %g,%g,%g", i, vec[0], vec[1], vec[2]);
900	#endif
901
902	double scale[3];
903	scale[0] = sqrt(vec[0]vec[0] + vec[1]vec[1] + vec[2]*vec[2]);
904	scale[1] = 1.0;
905	scale[2] = 1.0;
906
907	bondScales->InsertNextTupleValue(scale);
908	bondScales->InsertNextTupleValue(scale);
909
910	if (elements != NULL) {
911	int element = (int)elements->GetComponent(pts[0], 0);
912	#ifdef DEBUG
913	TRACE("Bond %d, elt 0: %d", i, element);
914	#endif
915	bondElements->InsertNextValue(element);
916	element = (int)elements->GetComponent(pts[1], 0);
917	#ifdef DEBUG
918	TRACE("Bond %d, elt 1: %d", i, element);
919	#endif
920	bondElements->InsertNextValue(element);
921	}
922	}
923	_bondPD->SetPoints(bondPoints);
924	if (elements != NULL)
925	_bondPD->GetPointData()->SetScalars(bondElements);
926	_bondPD->GetPointData()->AddArray(bondVectors);
927	_bondPD->GetPointData()->AddArray(bondScales);
928	}
929
930	void Molecule::addLabelArray(vtkDataSet *dataSet)
931	{
932	vtkDataArray *elements = NULL;
933	if (dataSet->GetPointData() != NULL &&
934	dataSet->GetPointData()->GetScalars() != NULL &&
935	strcmp(dataSet->GetPointData()->GetScalars()->GetName(), "element") == 0) {
936	elements = dataSet->GetPointData()->GetScalars();
937	} else {
938	WARN("Can't label atoms without an element array");
939	}
940
941	vtkSmartPointer<vtkStringArray> labelArray = vtkSmartPointer<vtkStringArray>::New();
942	labelArray->SetName("_atom_labels");
943	vtkPolyData *pd = vtkPolyData::SafeDownCast(dataSet);
944	if (pd == NULL) {
945	ERROR("DataSet not a PolyData");
946	return;
947	}
948	for (int i = 0; i < pd->GetNumberOfPoints(); i++) {
949	char buf[32];
950	if (elements != NULL) {
951	int elt = (int)elements->GetComponent(i, 0);
952	sprintf(buf, "%s%d", g_elementNames[elt], i);
953	} else {
954	sprintf(buf, "%d", i);
955	}
956	labelArray->InsertNextValue(buf);
957	}
958	dataSet->GetPointData()->AddArray(labelArray);
959	}
960
961	/**
962	* \brief Add a scalar array to dataSet with sizes for the elements
963	* specified in the "element" scalar array
964	*/
965	void Molecule::addRadiusArray(vtkDataSet *dataSet, AtomScaling scaling, double scaleFactor)
966	{
967	vtkDataArray *elements = NULL;
968	if (dataSet->GetPointData() != NULL &&
969	dataSet->GetPointData()->GetScalars() != NULL &&
970	strcmp(dataSet->GetPointData()->GetScalars()->GetName(), "element") == 0) {
971	elements = dataSet->GetPointData()->GetScalars();
972	} else if (scaling != NO_ATOM_SCALING) {
973	WARN("Can't use non-constant scaling without an element array");
974	}
975	const float *radiusSource = NULL;
976	switch (scaling) {
977	case VAN_DER_WAALS_RADIUS:
978	radiusSource = g_vdwRadii;
979	break;
980	case COVALENT_RADIUS:
981	radiusSource = g_covalentRadii;
982	break;
983	case ATOMIC_RADIUS:
984	radiusSource = g_atomicRadii;
985	break;
986	case NO_ATOM_SCALING:
987	default:
988	;
989	}
990	vtkPolyData *pd = vtkPolyData::SafeDownCast(dataSet);
991	if (pd == NULL) {
992	ERROR("DataSet not a PolyData");
993	return;
994	}
995	vtkDataArray *array = dataSet->GetPointData()->GetArray("_radii");
996	vtkSmartPointer<vtkFloatArray> radii;
997	if (array == NULL) {
998	radii = vtkSmartPointer<vtkFloatArray>::New();
999	radii->SetName("_radii");
1000	radii->SetNumberOfComponents(1);
1001	radii->SetNumberOfValues(pd->GetNumberOfPoints());
1002	} else {
1003	radii = vtkFloatArray::SafeDownCast(array);
1004	assert(radii != NULL);
1005	}
1006	for (int i = 0; i < pd->GetNumberOfPoints(); i++) {
1007	float value;
1008	if (elements != NULL && radiusSource != NULL) {
1009	int elt = (int)elements->GetComponent(i, 0);
1010	value = radiusSource[elt] * scaleFactor;
1011	} else {
1012	value = scaleFactor;
1013	}
1014	radii->SetValue(i, value);
1015	}
1016	radii->Modified();
1017	if (array == NULL) {
1018	dataSet->GetPointData()->AddArray(radii);
1019	}
1020	}
1021
1022	/**
1023	* \brief Create a color map to map atomic numbers to element colors
1024	*/
1025	ColorMap *Molecule::createElementColorMap()
1026	{
1027	ColorMap *elementCmap = new ColorMap("elementDefault");
1028	ColorMap::ControlPoint cp[NUM_ELEMENTS+1];
1029
1030	elementCmap->setNumberOfTableEntries(NUM_ELEMENTS+1);
1031	for (int i = 0; i <= NUM_ELEMENTS; i++) {
1032	cp[i].value = i/((double)NUM_ELEMENTS);
1033	for (int c = 0; c < 3; c++) {
1034	cp[i].color[c] = ((double)g_elementColors[i][c])/255.;
1035	}
1036	elementCmap->addControlPoint(cp[i]);
1037	}
1038	ColorMap::OpacityControlPoint ocp[2];
1039	ocp[0].value = 0;
1040	ocp[0].alpha = 1.0;
1041	ocp[1].value = 1.0;
1042	ocp[1].alpha = 1.0;
1043	elementCmap->addOpacityControlPoint(ocp[0]);
1044	elementCmap->addOpacityControlPoint(ocp[1]);
1045	elementCmap->build();
1046	double range[2];
1047	range[0] = 0;
1048	range[1] = NUM_ELEMENTS;
1049	elementCmap->getLookupTable()->SetRange(range);
1050
1051	return elementCmap;
1052	}
1053
1054	int Molecule::computeBonds(vtkDataSet *dataSet, float tolerance)
1055	{
1056	vtkPolyData *pd = vtkPolyData::SafeDownCast(dataSet);
1057	if (pd == NULL) {
1058	ERROR("DataSet not a PolyData");
1059	return -1;
1060	}
1061
1062	if (pd->GetNumberOfPoints() < 2 \|\|
1063	pd->GetNumberOfLines() > 0) {
1064	return 0;
1065	}
1066
1067	vtkDataArray *elements = NULL;
1068	if (pd->GetPointData() != NULL &&
1069	pd->GetPointData()->GetScalars() != NULL &&
1070	strcmp(pd->GetPointData()->GetScalars()->GetName(), "element") == 0) {
1071	elements = pd->GetPointData()->GetScalars();
1072	} else {
1073	TRACE("Can't compute bonds without an element array");
1074	return -1;
1075	}
1076	int numAtoms = pd->GetNumberOfPoints();
1077	vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
1078	for (int i = 0; i < numAtoms; i++) {
1079	double pt1[3];
1080	int elem1 = (int)elements->GetComponent(i, 0);
1081	float r1 = g_covalentRadii[elem1];
1082	pd->GetPoint(i, pt1);
1083	for (int j = i+1; j < numAtoms; j++) {
1084	float r2, dist;
1085	double pt2[3];
1086	int elem2 = (int)elements->GetComponent(j, 0);
1087
1088	if (elem1 == 1 && elem2 == 1)
1089	continue;
1090
1091	r2 = g_covalentRadii[elem2];
1092	pd->GetPoint(j, pt2);
1093
1094	double x = pt1[0] - pt2[0];
1095	double y = pt1[1] - pt2[1];
1096	double z = pt1[2] - pt2[2];
1097	dist = (float)sqrt(xx + yy + z*z);
1098	//TRACE("i=%d,j=%d elem1=%d,elem2=%d r1=%g,r2=%g, dist=%g", i, j, elem1, elem2, r1, r2, dist);
1099	if (dist < (r1 + r2 + tolerance)) {
1100	vtkIdType pts[2];
1101	pts[0] = i;
1102	pts[1] = j;
1103	cells->InsertNextCell(2, pts);
1104	}
1105	}
1106	}
1107
1108	if (cells->GetNumberOfCells() > 0) {
1109	pd->SetLines(cells);
1110	}
1111
1112	TRACE("Generated %d bonds", cells->GetNumberOfCells());
1113	return cells->GetNumberOfCells();
1114	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats: