Context Navigation

source: vtkvis/trunk/Molecule.cpp @ 5214

Last change on this file since 5214 was 5073, checked in by ldelgass, 9 years ago
merge r5072 from release branch
Property svn:eol-style set to `native`
File size: 34.5 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	#ifdef USE_VTK6
220	_bondMapper->SetInputData(_bondPD);
221	#else
222	_bondMapper->SetInputConnection(_bondPD->GetProducerPort());
223	#endif
224	_bondMapper->SetOrientationArray("bond_orientations");
225	_bondMapper->SetOrientationModeToDirection();
226	_bondMapper->OrientOn();
227	_bondMapper->SetScaleArray("bond_scales");
228	_bondMapper->SetScaleModeToScaleByVectorComponents();
229	_bondMapper->ScalingOn();
230	_bondMapper->ClampingOff();
231
232	_bondProp->SetMapper(_bondMapper);
233	getAssembly()->AddPart(_bondProp);
234	}
235	if (pd->GetNumberOfPoints() > 0) {
236	if (_labelHierarchy == NULL) {
237	_labelHierarchy = vtkSmartPointer<vtkPointSetToLabelHierarchy>::New();
238	}
239	if (_labelTransform == NULL) {
240	_labelTransform = vtkSmartPointer<vtkTransform>::New();
241	}
242	vtkSmartPointer<vtkTransformPolyDataFilter> transformFilter = vtkSmartPointer<vtkTransformPolyDataFilter>::New();
243	#ifdef USE_VTK6
244	transformFilter->SetInputData(pd);
245	#else
246	transformFilter->SetInput(pd);
247	#endif
248	transformFilter->SetTransform(_labelTransform);
249	_labelHierarchy->SetInputConnection(transformFilter->GetOutputPort());
250	_labelHierarchy->SetLabelArrayName("_atom_labels");
251	_labelHierarchy->GetTextProperty()->SetColor(0, 0, 0);
252	_labelMapper->SetInputConnection(_labelHierarchy->GetOutputPort());
253	_labelProp->SetMapper(_labelMapper);
254
255	// Atoms
256	if (_sphereSource == NULL) {
257	_sphereSource = vtkSmartPointer<vtkSphereSource>::New();
258	_sphereSource->SetRadius(1.0);
259	_sphereSource->SetThetaResolution(14);
260	_sphereSource->SetPhiResolution(14);
261	}
262
263	_atomMapper->SetSourceConnection(_sphereSource->GetOutputPort());
264	#ifdef USE_VTK6
265	_atomMapper->SetInputData(pd);
266	#else
267	_atomMapper->SetInputConnection(pd->GetProducerPort());
268	#endif
269	_atomMapper->SetScaleArray("_radii");
270	_atomMapper->SetScaleModeToScaleByMagnitude();
271	_atomMapper->ScalingOn();
272	_atomMapper->OrientOff();
273
274	_atomProp->SetMapper(_atomMapper);
275	getAssembly()->AddPart(_atomProp);
276	}
277	} else {
278	// DataSet is NOT a vtkPolyData
279	ERROR("DataSet is not a PolyData");
280	return;
281	}
282
283	setAtomLabelVisibility(_labelsOn);
284
285	if (pd->GetNumberOfPoints() > 0) {
286	_atomMapper->Update();
287	_labelMapper->Update();
288	}
289	if (pd->GetNumberOfLines() > 0) {
290	_bondMapper->Update();
291	}
292	}
293
294	void Molecule::updateLabelTransform()
295	{
296	if (_labelTransform != NULL && getAssembly() != NULL) {
297	_labelTransform->SetMatrix(getAssembly()->GetMatrix());
298	}
299	}
300
301	void Molecule::updateRanges(Renderer *renderer)
302	{
303	if (_dataSet == NULL) {
304	ERROR("called before setDataSet");
305	return;
306	}
307
308	if (renderer->getUseCumulativeRange()) {
309	renderer->getCumulativeDataRange(_dataRange,
310	_dataSet->getActiveScalarsName(),
311	1);
312	renderer->getCumulativeDataRange(_vectorMagnitudeRange,
313	_dataSet->getActiveVectorsName(),
314	3);
315	for (int i = 0; i < 3; i++) {
316	renderer->getCumulativeDataRange(_vectorComponentRange[i],
317	_dataSet->getActiveVectorsName(),
318	3, i);
319	}
320	} else {
321	_dataSet->getScalarRange(_dataRange);
322	_dataSet->getVectorRange(_vectorMagnitudeRange);
323	for (int i = 0; i < 3; i++) {
324	_dataSet->getVectorRange(_vectorComponentRange[i], i);
325	}
326	}
327
328	// Need to update color map ranges
329	double *rangePtr = _colorFieldRange;
330	if (_colorFieldRange[0] > _colorFieldRange[1]) {
331	rangePtr = NULL;
332	}
333	setColorMode(_colorMode, _colorFieldType, _colorFieldName.c_str(), rangePtr);
334	}
335
336	void Molecule::setColorMode(ColorMode mode)
337	{
338	_colorMode = mode;
339	if (_dataSet == NULL)
340	return;
341
342	switch (mode) {
343	case COLOR_BY_ELEMENTS: // Assume default scalar is "element" array
344	setColorMode(mode,
345	DataSet::POINT_DATA,
346	"element");
347	break;
348	case COLOR_BY_SCALAR:
349	setColorMode(mode,
350	_dataSet->getActiveScalarsType(),
351	_dataSet->getActiveScalarsName());
352	break;
353	case COLOR_BY_VECTOR_MAGNITUDE:
354	setColorMode(mode,
355	_dataSet->getActiveVectorsType(),
356	_dataSet->getActiveVectorsName());
357	break;
358	case COLOR_BY_VECTOR_X:
359	setColorMode(mode,
360	_dataSet->getActiveVectorsType(),
361	_dataSet->getActiveVectorsName());
362	break;
363	case COLOR_BY_VECTOR_Y:
364	setColorMode(mode,
365	_dataSet->getActiveVectorsType(),
366	_dataSet->getActiveVectorsName());
367	break;
368	case COLOR_BY_VECTOR_Z:
369	setColorMode(mode,
370	_dataSet->getActiveVectorsType(),
371	_dataSet->getActiveVectorsName());
372	break;
373	case COLOR_CONSTANT:
374	default:
375	setColorMode(mode, DataSet::POINT_DATA, NULL, NULL);
376	break;
377	}
378	}
379
380	void Molecule::setColorMode(ColorMode mode,
381	const char *name, double range[2])
382	{
383	if (_dataSet == NULL)
384	return;
385	DataSet::DataAttributeType type = DataSet::POINT_DATA;
386	int numComponents = 1;
387	if (name != NULL && strlen(name) > 0 &&
388	!_dataSet->getFieldInfo(name, &type, &numComponents)) {
389	ERROR("Field not found: %s", name);
390	return;
391	}
392	setColorMode(mode, type, name, range);
393	}
394
395	void Molecule::setColorMode(ColorMode mode, DataSet::DataAttributeType type,
396	const char *name, double range[2])
397	{
398	_colorMode = mode;
399	_colorFieldType = type;
400	if (name == NULL)
401	_colorFieldName.clear();
402	else
403	_colorFieldName = name;
404	if (range == NULL) {
405	_colorFieldRange[0] = DBL_MAX;
406	_colorFieldRange[1] = -DBL_MAX;
407	} else {
408	memcpy(_colorFieldRange, range, sizeof(double)*2);
409	}
410
411	if (_dataSet == NULL \|\| (_atomMapper == NULL && _bondMapper == NULL))
412	return;
413
414	switch (type) {
415	case DataSet::POINT_DATA:
416	if (_atomMapper != NULL)
417	_atomMapper->SetScalarModeToUsePointFieldData();
418	if (_bondMapper != NULL)
419	_bondMapper->SetScalarModeToUsePointFieldData();
420	break;
421	case DataSet::CELL_DATA:
422	if (_atomMapper != NULL)
423	_atomMapper->SetScalarModeToUseCellFieldData();
424	if (_bondMapper != NULL)
425	_bondMapper->SetScalarModeToUseCellFieldData();
426	break;
427	default:
428	ERROR("Unsupported DataAttributeType: %d", type);
429	return;
430	}
431
432	if (name != NULL && strlen(name) > 0) {
433	if (_atomMapper != NULL)
434	_atomMapper->SelectColorArray(name);
435	if (_bondMapper != NULL)
436	_bondMapper->SelectColorArray(name);
437	} else {
438	if (_atomMapper != NULL)
439	_atomMapper->SetScalarModeToDefault();
440	if (_bondMapper != NULL)
441	_bondMapper->SetScalarModeToDefault();
442	}
443
444	if (_lut != NULL) {
445	if (range != NULL) {
446	_lut->SetRange(range);
447	} else if (mode == COLOR_BY_ELEMENTS) {
448	double range[2];
449	range[0] = 0;
450	range[1] = NUM_ELEMENTS;
451	_lut->SetRange(range);
452	} else if (name != NULL && strlen(name) > 0) {
453	double r[2];
454	int comp = -1;
455	if (mode == COLOR_BY_VECTOR_X)
456	comp = 0;
457	else if (mode == COLOR_BY_VECTOR_Y)
458	comp = 1;
459	else if (mode == COLOR_BY_VECTOR_Z)
460	comp = 2;
461
462	if (_renderer->getUseCumulativeRange()) {
463	int numComponents;
464	if (!_dataSet->getFieldInfo(name, type, &numComponents)) {
465	ERROR("Field not found: %s, type: %d", name, type);
466	return;
467	} else if (numComponents < comp+1) {
468	ERROR("Request for component %d in field with %d components",
469	comp, numComponents);
470	return;
471	}
472	_renderer->getCumulativeDataRange(r, name, type, numComponents, comp);
473	} else {
474	_dataSet->getDataRange(r, name, type, comp);
475	}
476	_lut->SetRange(r);
477	}
478	}
479
480	switch (mode) {
481	case COLOR_BY_ELEMENTS:
482	if (_atomMapper != NULL)
483	_atomMapper->ScalarVisibilityOn();
484	if (_bondMapper != NULL)
485	_bondMapper->ScalarVisibilityOn();
486	break;
487	case COLOR_BY_SCALAR:
488	if (_atomMapper != NULL)
489	_atomMapper->ScalarVisibilityOn();
490	if (_bondMapper != NULL)
491	_bondMapper->ScalarVisibilityOn();
492	break;
493	case COLOR_BY_VECTOR_MAGNITUDE:
494	if (_atomMapper != NULL)
495	_atomMapper->ScalarVisibilityOn();
496	if (_bondMapper != NULL)
497	_bondMapper->ScalarVisibilityOn();
498	if (_lut != NULL) {
499	_lut->SetVectorModeToMagnitude();
500	}
501	break;
502	case COLOR_BY_VECTOR_X:
503	if (_atomMapper != NULL)
504	_atomMapper->ScalarVisibilityOn();
505	if (_bondMapper != NULL)
506	_bondMapper->ScalarVisibilityOn();
507	if (_lut != NULL) {
508	_lut->SetVectorModeToComponent();
509	_lut->SetVectorComponent(0);
510	}
511	break;
512	case COLOR_BY_VECTOR_Y:
513	if (_atomMapper != NULL)
514	_atomMapper->ScalarVisibilityOn();
515	if (_bondMapper != NULL)
516	_bondMapper->ScalarVisibilityOn();
517	if (_lut != NULL) {
518	_lut->SetVectorModeToComponent();
519	_lut->SetVectorComponent(1);
520	}
521	break;
522	case COLOR_BY_VECTOR_Z:
523	if (_atomMapper != NULL)
524	_atomMapper->ScalarVisibilityOn();
525	if (_bondMapper != NULL)
526	_bondMapper->ScalarVisibilityOn();
527	if (_lut != NULL) {
528	_lut->SetVectorModeToComponent();
529	_lut->SetVectorComponent(2);
530	}
531	break;
532	case COLOR_CONSTANT:
533	default:
534	if (_atomMapper != NULL)
535	_atomMapper->ScalarVisibilityOff();
536	if (_bondMapper != NULL)
537	_bondMapper->ScalarVisibilityOff();
538	break;
539	}
540	}
541
542	/**
543	* \brief Called when the color map has been edited
544	*/
545	void Molecule::updateColorMap()
546	{
547	setColorMap(_colorMap);
548	}
549
550	/**
551	* \brief Associate a colormap lookup table with the DataSet
552	*/
553	void Molecule::setColorMap(ColorMap *cmap)
554	{
555	if (cmap == NULL)
556	return;
557
558	_colorMap = cmap;
559
560	if (_lut == NULL) {
561	_lut = vtkSmartPointer<vtkLookupTable>::New();
562	if (_atomMapper != NULL) {
563	_atomMapper->UseLookupTableScalarRangeOn();
564	_atomMapper->SetLookupTable(_lut);
565	}
566	if (_bondMapper != NULL) {
567	_bondMapper->UseLookupTableScalarRangeOn();
568	_bondMapper->SetLookupTable(_lut);
569	}
570	_lut->DeepCopy(cmap->getLookupTable());
571	switch (_colorMode) {
572	case COLOR_BY_ELEMENTS: {
573	double range[2];
574	range[0] = 0;
575	range[1] = NUM_ELEMENTS;
576	_lut->SetRange(range);
577	}
578	break;
579	case COLOR_CONSTANT:
580	case COLOR_BY_SCALAR:
581	_lut->SetRange(_dataRange);
582	break;
583	case COLOR_BY_VECTOR_MAGNITUDE:
584	_lut->SetRange(_vectorMagnitudeRange);
585	break;
586	case COLOR_BY_VECTOR_X:
587	_lut->SetRange(_vectorComponentRange[0]);
588	break;
589	case COLOR_BY_VECTOR_Y:
590	_lut->SetRange(_vectorComponentRange[1]);
591	break;
592	case COLOR_BY_VECTOR_Z:
593	_lut->SetRange(_vectorComponentRange[2]);
594	break;
595	default:
596	break;
597	}
598	} else {
599	double range[2];
600	_lut->GetTableRange(range);
601	_lut->DeepCopy(cmap->getLookupTable());
602	_lut->SetRange(range);
603	_lut->Modified();
604	}
605	}
606
607	void Molecule::setAtomLabelField(const char *fieldName)
608	{
609	if (_labelHierarchy != NULL) {
610	if (strcmp(fieldName, "default") == 0) {
611	_labelHierarchy->SetLabelArrayName("_atom_labels");
612	} else {
613	_labelHierarchy->SetLabelArrayName(fieldName);
614	}
615	}
616	}
617
618	/**
619	* \brief Turn on/off rendering of atom labels
620	*/
621	void Molecule::setAtomLabelVisibility(bool state)
622	{
623	_labelsOn = state;
624	if (_labelProp != NULL) {
625	_labelProp->SetVisibility((state ? 1 : 0));
626	}
627	}
628
629	/**
630	* \brief Turn on/off rendering of the atoms
631	*/
632	void Molecule::setAtomVisibility(bool state)
633	{
634	if (_atomProp != NULL) {
635	_atomProp->SetVisibility((state ? 1 : 0));
636	}
637	}
638
639	/**
640	* \brief Turn on/off rendering of the bonds
641	*/
642	void Molecule::setBondVisibility(bool state)
643	{
644	if (_bondProp != NULL) {
645	_bondProp->SetVisibility((state ? 1 : 0));
646	}
647	}
648
649	/**
650	* \brief Toggle visibility of the prop
651	*/
652	void Molecule::setVisibility(bool state)
653	{
654	GraphicsObject::setVisibility(state);
655	if (_labelProp != NULL) {
656	if (!state)
657	_labelProp->SetVisibility(0);
658	else
659	setAtomLabelVisibility(_labelsOn);
660	}
661	}
662
663	/**
664	* \brief Set opacity of molecule
665	*/
666	void Molecule::setOpacity(double opacity)
667	{
668	GraphicsObject::setOpacity(opacity);
669	if (_labelMapper != NULL) {
670	_labelMapper->SetBackgroundOpacity(opacity);
671	}
672	}
673
674	void Molecule::setAtomQuality(double quality)
675	{
676	if (_sphereSource == NULL)
677	return;
678
679	if (quality > 10.0)
680	quality = 10.0;
681
682	int thetaRes = (int)(quality * 14.0);
683	int phiRes = (int)(quality * 14.0);
684	if (thetaRes < 4) thetaRes = 4;
685	if (phiRes < 3) phiRes = 3;
686
687	_sphereSource->SetThetaResolution(thetaRes);
688	_sphereSource->SetPhiResolution(phiRes);
689
690	if (_atomMapper != NULL) {
691	_atomMapper->Modified();
692	_atomMapper->Update();
693	}
694	}
695
696	void Molecule::setBondQuality(double quality)
697	{
698	if (_cylinderSource == NULL)
699	return;
700
701	if (quality > 10.0)
702	quality = 10.0;
703
704	int res = (int)(quality * 12.0);
705	if (res < 3) res = 3;
706
707	_cylinderSource->SetResolution(res);
708
709	if (_bondMapper != NULL) {
710	_bondMapper->Modified();
711	_bondMapper->Update();
712	}
713	}
714
715	void Molecule::setBondStyle(BondStyle style)
716	{
717	switch (style) {
718	case BOND_STYLE_CYLINDER:
719	if (_bondProp != NULL) {
720	_bondProp->GetProperty()->SetLineWidth(_edgeWidth);
721	_bondProp->GetProperty()->SetLighting(_lighting ? 1 : 0);
722	}
723	if (_bondMapper != NULL && _cylinderTrans != NULL &&
724	_bondMapper->GetInputConnection(1, 0) != _cylinderTrans->GetOutputPort()) {
725	_cylinderTrans->Modified();
726	_bondMapper->SetSourceConnection(_cylinderTrans->GetOutputPort());
727	_bondMapper->Modified();
728	}
729	break;
730	case BOND_STYLE_LINE:
731	if (_bondProp != NULL) {
732	_bondProp->GetProperty()->LightingOff();
733	}
734	if (_bondMapper != NULL && _lineSource != NULL &&
735	_bondMapper->GetInputConnection(1, 0) != _lineSource->GetOutputPort()) {
736	_lineSource->Modified();
737	_bondMapper->SetSourceConnection(_lineSource->GetOutputPort());
738	_bondMapper->Modified();
739	}
740	break;
741	default:
742	WARN("Unknown bond style");
743	}
744	}
745
746	/**
747	* \brief Set constant bond color
748	*/
749	void Molecule::setBondColor(float color[3])
750	{
751	_bondColor[0] = color[0];
752	_bondColor[1] = color[1];
753	_bondColor[2] = color[2];
754	if (_bondProp != NULL) {
755	_bondProp->GetProperty()->SetColor(_bondColor[0], _bondColor[1], _bondColor[2]);
756	}
757	}
758
759	/**
760	* \brief Set mode to determine how bonds are colored
761	*/
762	void Molecule::setBondColorMode(BondColorMode mode)
763	{
764	if (_bondMapper == NULL) return;
765
766	switch (mode) {
767	case BOND_COLOR_BY_ELEMENTS:
768	_bondMapper->ScalarVisibilityOn();
769	break;
770	case BOND_COLOR_CONSTANT:
771	_bondMapper->ScalarVisibilityOff();
772	break;
773	default:
774	WARN("Unknown bond color mode");
775	}
776	}
777
778	/**
779	* \brief Set a group of world coordinate planes to clip rendering
780	*
781	* Passing NULL for planes will remove all cliping planes
782	*/
783	void Molecule::setClippingPlanes(vtkPlaneCollection *planes)
784	{
785	if (_atomMapper != NULL) {
786	_atomMapper->SetClippingPlanes(planes);
787	}
788	if (_bondMapper != NULL) {
789	_bondMapper->SetClippingPlanes(planes);
790	}
791	}
792
793	/**
794	* \brief Set the radius type used for scaling atoms
795	*/
796	void Molecule::setAtomScaling(AtomScaling state)
797	{
798	if (state == _atomScaling)
799	return;
800
801	_atomScaling = state;
802	if (_dataSet != NULL) {
803	vtkDataSet *ds = _dataSet->getVtkDataSet();
804	addRadiusArray(ds, _atomScaling, _radiusScale);
805	if (_atomMapper != NULL) {
806	assert(ds->GetPointData() != NULL &&
807	ds->GetPointData()->GetArray("_radii") != NULL);
808	_atomMapper->SetScaleModeToScaleByMagnitude();
809	_atomMapper->SetScaleArray("_radii");
810	_atomMapper->ScalingOn();
811	}
812	}
813	}
814
815	/**
816	* \brief Set the constant radius scaling factor for atoms. This
817	* can be used to convert from Angstroms to atom coordinates units.
818	*/
819	void Molecule::setAtomRadiusScale(double scale)
820	{
821	if (scale == _radiusScale)
822	return;
823
824	_radiusScale = scale;
825	if (_dataSet != NULL) {
826	vtkDataSet *ds = _dataSet->getVtkDataSet();
827	addRadiusArray(ds, _atomScaling, _radiusScale);
828	}
829	}
830
831	/**
832	* \brief Set the constant radius scaling factor for bonds.
833	*/
834	void Molecule::setBondRadiusScale(double scale)
835	{
836	if (_cylinderSource != NULL &&
837	_cylinderSource->GetRadius() != scale) {
838	_cylinderSource->SetRadius(scale);
839	// Workaround bug with source modification not causing
840	// mapper to be updated
841	if (_bondMapper != NULL) {
842	_bondMapper->Modified();
843	}
844	}
845	}
846
847	void Molecule::setupBondPolyData()
848	{
849	if (_dataSet == NULL)
850	return;
851	if (_bondPD == NULL) {
852	_bondPD = vtkSmartPointer<vtkPolyData>::New();
853	} else {
854	_bondPD->Initialize();
855	}
856	vtkPolyData *pd = vtkPolyData::SafeDownCast(_dataSet->getVtkDataSet());
857	if (pd == NULL)
858	return;
859	vtkCellArray *lines = pd->GetLines();
860	lines->InitTraversal();
861	vtkIdType npts, *pts;
862	vtkSmartPointer<vtkPoints> bondPoints = vtkSmartPointer<vtkPoints>::New();
863	vtkSmartPointer<vtkIntArray> bondElements = vtkSmartPointer<vtkIntArray>::New();
864	vtkSmartPointer<vtkDoubleArray> bondVectors = vtkSmartPointer<vtkDoubleArray>::New();
865	vtkSmartPointer<vtkDoubleArray> bondScales = vtkSmartPointer<vtkDoubleArray>::New();
866	bondElements->SetName("element");
867	bondVectors->SetName("bond_orientations");
868	bondVectors->SetNumberOfComponents(3);
869	bondScales->SetName("bond_scales");
870	bondScales->SetNumberOfComponents(3);
871	vtkDataArray *elements = NULL;
872	if (pd->GetPointData() != NULL &&
873	pd->GetPointData()->GetScalars() != NULL &&
874	strcmp(pd->GetPointData()->GetScalars()->GetName(), "element") == 0) {
875	elements = pd->GetPointData()->GetScalars();
876	}
877	for (int i = 0; lines->GetNextCell(npts, pts); i++) {
878	assert(npts == 2);
879	double pt0[3], pt1[3];
880	double newPt0[3], newPt1[3];
881	//double *pt0 = pd->GetPoint(pts[0]);
882	//double *pt1 = pd->GetPoint(pts[1]);
883	pd->GetPoint(pts[0], pt0);
884	pd->GetPoint(pts[1], pt1);
885	double center[3];
886
887	for (int j = 0; j < 3; j++)
888	center[j] = pt0[j] + (pt1[j] - pt0[j]) * 0.5;
889	for (int j = 0; j < 3; j++)
890	newPt0[j] = pt0[j] + (pt1[j] - pt0[j]) * 0.25;
891	for (int j = 0; j < 3; j++)
892	newPt1[j] = pt0[j] + (pt1[j] - pt0[j]) * 0.75;
893
894	bondPoints->InsertNextPoint(newPt0);
895	bondPoints->InsertNextPoint(newPt1);
896
897	#ifdef DEBUG
898	TRACE("Bond %d: (%g,%g,%g)-(%g,%g,%g)-(%g,%g,%g)", i,
899	pt0[0], pt0[1], pt0[2],
900	center[0], center[1], center[2],
901	pt1[0], pt1[1], pt1[2]);
902	#endif
903
904	double vec[3];
905	for (int j = 0; j < 3; j++)
906	vec[j] = center[j] - pt0[j];
907
908	bondVectors->InsertNextTupleValue(vec);
909	bondVectors->InsertNextTupleValue(vec);
910	#ifdef DEBUG
911	TRACE("Bond %d, vec: %g,%g,%g", i, vec[0], vec[1], vec[2]);
912	#endif
913
914	double scale[3];
915	scale[0] = sqrt(vec[0]vec[0] + vec[1]vec[1] + vec[2]*vec[2]);
916	scale[1] = 1.0;
917	scale[2] = 1.0;
918
919	bondScales->InsertNextTupleValue(scale);
920	bondScales->InsertNextTupleValue(scale);
921
922	if (elements != NULL) {
923	int element = (int)elements->GetComponent(pts[0], 0);
924	#ifdef DEBUG
925	TRACE("Bond %d, elt 0: %d", i, element);
926	#endif
927	bondElements->InsertNextValue(element);
928	element = (int)elements->GetComponent(pts[1], 0);
929	#ifdef DEBUG
930	TRACE("Bond %d, elt 1: %d", i, element);
931	#endif
932	bondElements->InsertNextValue(element);
933	}
934	}
935	_bondPD->SetPoints(bondPoints);
936	if (elements != NULL)
937	_bondPD->GetPointData()->SetScalars(bondElements);
938	_bondPD->GetPointData()->AddArray(bondVectors);
939	_bondPD->GetPointData()->AddArray(bondScales);
940	}
941
942	void Molecule::addLabelArray(vtkDataSet *dataSet)
943	{
944	vtkDataArray *elements = NULL;
945	if (dataSet->GetPointData() != NULL &&
946	dataSet->GetPointData()->GetScalars() != NULL &&
947	strcmp(dataSet->GetPointData()->GetScalars()->GetName(), "element") == 0) {
948	elements = dataSet->GetPointData()->GetScalars();
949	} else {
950	WARN("Can't label atoms without an element array");
951	}
952
953	vtkSmartPointer<vtkStringArray> labelArray = vtkSmartPointer<vtkStringArray>::New();
954	labelArray->SetName("_atom_labels");
955	vtkPolyData *pd = vtkPolyData::SafeDownCast(dataSet);
956	if (pd == NULL) {
957	ERROR("DataSet not a PolyData");
958	return;
959	}
960	for (int i = 0; i < pd->GetNumberOfPoints(); i++) {
961	char buf[32];
962	if (elements != NULL) {
963	int elt = (int)elements->GetComponent(i, 0);
964	sprintf(buf, "%s%d", g_elementNames[elt], i);
965	} else {
966	sprintf(buf, "%d", i);
967	}
968	labelArray->InsertNextValue(buf);
969	}
970	dataSet->GetPointData()->AddArray(labelArray);
971	}
972
973	/**
974	* \brief Add a scalar array to dataSet with sizes for the elements
975	* specified in the "element" scalar array
976	*/
977	void Molecule::addRadiusArray(vtkDataSet *dataSet, AtomScaling scaling, double scaleFactor)
978	{
979	vtkDataArray *elements = NULL;
980	if (dataSet->GetPointData() != NULL &&
981	dataSet->GetPointData()->GetScalars() != NULL &&
982	strcmp(dataSet->GetPointData()->GetScalars()->GetName(), "element") == 0) {
983	elements = dataSet->GetPointData()->GetScalars();
984	} else if (scaling != NO_ATOM_SCALING) {
985	WARN("Can't use non-constant scaling without an element array");
986	}
987	const float *radiusSource = NULL;
988	switch (scaling) {
989	case VAN_DER_WAALS_RADIUS:
990	radiusSource = g_vdwRadii;
991	break;
992	case COVALENT_RADIUS:
993	radiusSource = g_covalentRadii;
994	break;
995	case ATOMIC_RADIUS:
996	radiusSource = g_atomicRadii;
997	break;
998	case NO_ATOM_SCALING:
999	default:
1000	;
1001	}
1002	vtkPolyData *pd = vtkPolyData::SafeDownCast(dataSet);
1003	if (pd == NULL) {
1004	ERROR("DataSet not a PolyData");
1005	return;
1006	}
1007	vtkDataArray *array = dataSet->GetPointData()->GetArray("_radii");
1008	vtkSmartPointer<vtkFloatArray> radii;
1009	if (array == NULL) {
1010	radii = vtkSmartPointer<vtkFloatArray>::New();
1011	radii->SetName("_radii");
1012	radii->SetNumberOfComponents(1);
1013	radii->SetNumberOfValues(pd->GetNumberOfPoints());
1014	} else {
1015	radii = vtkFloatArray::SafeDownCast(array);
1016	assert(radii != NULL);
1017	}
1018	for (int i = 0; i < pd->GetNumberOfPoints(); i++) {
1019	float value;
1020	if (elements != NULL && radiusSource != NULL) {
1021	int elt = (int)elements->GetComponent(i, 0);
1022	value = radiusSource[elt] * scaleFactor;
1023	} else {
1024	value = scaleFactor;
1025	}
1026	radii->SetValue(i, value);
1027	}
1028	radii->Modified();
1029	if (array == NULL) {
1030	dataSet->GetPointData()->AddArray(radii);
1031	}
1032	}
1033
1034	/**
1035	* \brief Create a color map to map atomic numbers to element colors
1036	*/
1037	ColorMap *Molecule::createElementColorMap()
1038	{
1039	ColorMap *elementCmap = new ColorMap("elementDefault");
1040	ColorMap::ControlPoint cp[NUM_ELEMENTS+1];
1041
1042	elementCmap->setNumberOfTableEntries(NUM_ELEMENTS+1);
1043	for (int i = 0; i <= NUM_ELEMENTS; i++) {
1044	cp[i].value = i/((double)NUM_ELEMENTS);
1045	for (int c = 0; c < 3; c++) {
1046	cp[i].color[c] = ((double)g_elementColors[i][c])/255.;
1047	}
1048	elementCmap->addControlPoint(cp[i]);
1049	}
1050	ColorMap::OpacityControlPoint ocp[2];
1051	ocp[0].value = 0;
1052	ocp[0].alpha = 1.0;
1053	ocp[1].value = 1.0;
1054	ocp[1].alpha = 1.0;
1055	elementCmap->addOpacityControlPoint(ocp[0]);
1056	elementCmap->addOpacityControlPoint(ocp[1]);
1057	elementCmap->build();
1058	double range[2];
1059	range[0] = 0;
1060	range[1] = NUM_ELEMENTS;
1061	elementCmap->getLookupTable()->SetRange(range);
1062
1063	return elementCmap;
1064	}
1065
1066	int Molecule::computeBonds(vtkDataSet *dataSet, float tolerance)
1067	{
1068	vtkPolyData *pd = vtkPolyData::SafeDownCast(dataSet);
1069	if (pd == NULL) {
1070	ERROR("DataSet not a PolyData");
1071	return -1;
1072	}
1073
1074	if (pd->GetNumberOfPoints() < 2 \|\|
1075	pd->GetNumberOfLines() > 0) {
1076	return 0;
1077	}
1078
1079	vtkDataArray *elements = NULL;
1080	if (pd->GetPointData() != NULL &&
1081	pd->GetPointData()->GetScalars() != NULL &&
1082	strcmp(pd->GetPointData()->GetScalars()->GetName(), "element") == 0) {
1083	elements = pd->GetPointData()->GetScalars();
1084	} else {
1085	TRACE("Can't compute bonds without an element array");
1086	return -1;
1087	}
1088	int numAtoms = pd->GetNumberOfPoints();
1089	vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
1090	for (int i = 0; i < numAtoms; i++) {
1091	double pt1[3];
1092	int elem1 = (int)elements->GetComponent(i, 0);
1093	float r1 = g_covalentRadii[elem1];
1094	pd->GetPoint(i, pt1);
1095	for (int j = i+1; j < numAtoms; j++) {
1096	float r2, dist;
1097	double pt2[3];
1098	int elem2 = (int)elements->GetComponent(j, 0);
1099
1100	if (elem1 == 1 && elem2 == 1)
1101	continue;
1102
1103	r2 = g_covalentRadii[elem2];
1104	pd->GetPoint(j, pt2);
1105
1106	double x = pt1[0] - pt2[0];
1107	double y = pt1[1] - pt2[1];
1108	double z = pt1[2] - pt2[2];
1109	dist = (float)sqrt(xx + yy + z*z);
1110	//TRACE("i=%d,j=%d elem1=%d,elem2=%d r1=%g,r2=%g, dist=%g", i, j, elem1, elem2, r1, r2, dist);
1111	if (dist < (r1 + r2 + tolerance)) {
1112	vtkIdType pts[2];
1113	pts[0] = i;
1114	pts[1] = j;
1115	cells->InsertNextCell(2, pts);
1116	}
1117	}
1118	}
1119
1120	if (cells->GetNumberOfCells() > 0) {
1121	pd->SetLines(cells);
1122	}
1123
1124	TRACE("Generated %d bonds", cells->GetNumberOfCells());
1125	return cells->GetNumberOfCells();
1126	}

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

Download in other formats: