Context Navigation

source: vtkvis/trunk/Molecule.cpp @ 4608

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

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

Download in other formats: