Changeset 2349 for trunk/packages/vizservers

Timestamp:

Aug 12, 2011 1:05:54 PM (13 years ago)

Author:

ldelgass

Message:

Add disk and filled regular polygon streamline seed sources, allow rotating
polygon seed about its normal. Also, add implementation of relative quat
rotations of camera and add flag to control relative v. absolute.

Location:

trunk/packages/vizservers/vtkvis

Files:

• RpStreamlines.cpp (modified) (12 diffs)
• RpStreamlines.h (modified) (2 diffs)
• RpVtkRenderer.cpp (modified) (9 diffs)
• RpVtkRenderer.h (modified) (4 diffs)
• RpVtkRendererCmd.cpp (modified) (4 diffs)
• protocol.txt (modified) (2 diffs)

trunk/packages/vizservers/vtkvis/RpStreamlines.cpp

@@ -9 +9 @@
 #include <ctime>
 #include <cfloat>
+#include <cmath>
+
+#include <vtkMath.h>
 #include <vtkActor.h>
 #include <vtkProperty.h>
@@ -21 +24 @@
 #include <vtkTubeFilter.h>
 #include <vtkRibbonFilter.h>
+#include <vtkTransform.h>
+#include <vtkTransformPolyDataFilter.h>
 
 #include "RpStreamlines.h"
@@ -36 +41 @@
     _seedColor[1] = 1.0f;
     _seedColor[2] = 1.0f;
+    vtkMath::RandomSeed((int)time(NULL));
+    srand((unsigned int)time(NULL));
 }
 
@@ -92 +99 @@
 }
 
+/**
+ * \brief Get a pseudo-random number in range [min,max]
+ */
+double Streamlines::getRandomNum(double min, double max)
+{
+#if 1
+    return vtkMath::Random(min, max);
+#else
+    int r = rand();
+    return (min + ((double)r / RAND_MAX) * (max - min));
+#endif
+}
+
+/**
+ * \brief Get a random 3D point within an AABB
+ *
+ * \param[out] pt The random point
+ * \param[in] bounds The bounds of the AABB
+ */
 void Streamlines::getRandomPoint(double pt[3], const double bounds[6])
 {
-    int r = rand();
-    pt[0] = bounds[0] + ((double)r / RAND_MAX) * (bounds[1] - bounds[0]);
-    r = rand();
-    pt[1] = bounds[2] + ((double)r / RAND_MAX) * (bounds[3] - bounds[2]);
-    r = rand();
-    pt[2] = bounds[4] + ((double)r / RAND_MAX) * (bounds[5] - bounds[4]);
-}
-
-void Streamlines::getRandomCellPt(vtkDataSet *ds, double pt[3])
+    pt[0] = getRandomNum(bounds[0], bounds[1]);
+    pt[1] = getRandomNum(bounds[2], bounds[3]);
+    pt[2] = getRandomNum(bounds[4], bounds[5]);
+}
+
+/**
+ * \brief Get a random point within a triangle (including edges)
+ *
+ * \param[out] pt The random point
+ * \param[in] v1 Triangle vertex 1
+ * \param[in] v2 Triangle vertex 2
+ * \param[in] v3 Triangle vertex 3
+ */
+void Streamlines::getRandomPointInTriangle(double pt[3],
+                                           const double v1[3],
+                                           const double v2[3],
+                                           const double v3[3])
+{
+    // Choose random barycentric coordinates
+    double bary[3];
+    bary[0] = getRandomNum(0, 1);
+    bary[1] = getRandomNum(0, 1);
+    if (bary[0] + bary[1] > 1.0) {
+        bary[0] = 1.0 - bary[0];
+        bary[1] = 1.0 - bary[1];
+    }
+    bary[2] = 1.0 - bary[0] - bary[1];
+
+    TRACE("bary %g %g %g", bary[0], bary[1], bary[2]);
+    // Convert to cartesian coords
+    for (int i = 0; i < 3; i++) {
+        pt[i] = v1[i] * bary[0] + v2[i] * bary[1] + v3[i] * bary[2];
+    }
+}
+
+/**
+ * \brief Get a random point on a line segment (including endpoints)
+ */
+void Streamlines::getRandomPointOnLineSegment(double pt[3],
+                                              const double endpt[3],
+                                              const double endpt2[3])
+{
+    double ratio = getRandomNum(0, 1);
+    pt[0] = endpt[0] + ratio * (endpt2[0] - endpt[0]);
+    pt[1] = endpt[1] + ratio * (endpt2[1] - endpt[1]);
+    pt[2] = endpt[2] + ratio * (endpt2[2] - endpt[2]);
+}
+
+/**
+ * \brief Get a random point within a vtkDataSet's mesh
+ *
+ * Note: This currently doesn't give a uniform distribution of
+ * points in space and can generate points outside the mesh
+ */
+void Streamlines::getRandomCellPt(double pt[3], vtkDataSet *ds)
 {
     int numCells = (int)ds->GetNumberOfCells();
+    // XXX: Not uniform distribution (shouldn't use mod, and assumes
+    // all cells are equal area/volume)
     int cell = rand() % numCells;
     double bounds[6];
     ds->GetCellBounds(cell, bounds);
-    int r = rand();
-    pt[0] = bounds[0] + ((double)r / RAND_MAX) * (bounds[1] - bounds[0]);
-    r = rand();
-    pt[1] = bounds[2] + ((double)r / RAND_MAX) * (bounds[3] - bounds[2]);
-    r = rand();
-    pt[2] = bounds[4] + ((double)r / RAND_MAX) * (bounds[5] - bounds[4]);
+    // Note: point is inside AABB of cell, but may be outside the cell
+    getRandomPoint(pt, bounds);
 }
 
@@ -164 +234 @@
 
     _streamTracer->SetInput(ds);
-
-    // Set up seed source object
-    vtkSmartPointer<vtkPolyData> seed = vtkSmartPointer<vtkPolyData>::New();
-    vtkSmartPointer<vtkPoints> pts = vtkSmartPointer<vtkPoints>::New();
-    vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
-
-    int numPoints = 200;
-    srand((unsigned int)time(NULL));
-    for (int i = 0; i < numPoints; i++) {
-        double pt[3];
-        getRandomCellPt(ds, pt);
-        TRACE("Seed pt: %g %g %g", pt[0], pt[1], pt[2]);
-        pts->InsertNextPoint(pt);
-        cells->InsertNextCell(1);
-        cells->InsertCellPoint(i);
-    }
-
-    seed->SetPoints(pts);
-    seed->SetVerts(cells);
-
-    TRACE("Seed points: %d", seed->GetNumberOfPoints());
-
     _streamTracer->SetMaximumPropagation(maxBound);
-    _streamTracer->SetSource(seed);
 
     if (_pdMapper == NULL) {
@@ -193 +240 @@
         _pdMapper->SetResolveCoincidentTopologyToPolygonOffset();
         _pdMapper->ScalarVisibilityOn();
-        // _pdMapper->ScalarVisibilityOff();
     }
     if (_seedMapper == NULL) {
         _seedMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
         _seedMapper->SetResolveCoincidentTopologyToPolygonOffset();
-        _seedMapper->SetInput(seed);
-    }
+    }
+
+    // Set up seed source object
+    setSeedToRandomPoints(200);
 
     switch (_lineType) {
@@ -262 +310 @@
  * \brief Use randomly distributed seed points
  *
+ * Note: The current implementation doesn't give a uniform 
+ * distribution of points, and points outside the mesh bounds
+ * may be generated
+ *
  * \param[in] numPoints Number of random seed points to generate
  */
@@ -272 +324 @@
         vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
 
-        srand((unsigned int)time(NULL));
         for (int i = 0; i < numPoints; i++) {
             double pt[3];
-            getRandomCellPt(_dataSet->getVtkDataSet(), pt);
+            getRandomCellPt(pt, _dataSet->getVtkDataSet());
             TRACE("Seed pt: %g %g %g", pt[0], pt[1], pt[2]);
             pts->InsertNextPoint(pt);
@@ -354 +405 @@
 
 /**
+ * \brief Create seed points inside a disk with an optional hole
+ *
+ * \param[in] center Center point of disk
+ * \param[in] normal Normal vector to orient disk
+ * \param[in] radius Radius of disk
+ * \param[in] innerRadius Radius of hole at center of disk
+ * \param[in] numPoints Number of random points to generate
+ */
+void Streamlines::setSeedToDisk(double center[3],
+                                double normal[3],
+                                double radius,
+                                double innerRadius,
+                                int numPoints)
+{
+    if (_streamTracer != NULL) {
+        // Set up seed source object
+        vtkSmartPointer<vtkPolyData> seed = vtkSmartPointer<vtkPolyData>::New();
+        vtkSmartPointer<vtkPoints> pts = vtkSmartPointer<vtkPoints>::New();
+        vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
+
+        // The following code is based on vtkRegularPolygonSource::RequestData
+
+        double px[3];
+        double py[3];
+        double axis[3] = {1., 0., 0.};
+
+        if (vtkMath::Normalize(normal) == 0.0) {
+            normal[0] = 0.0;
+            normal[1] = 0.0;
+            normal[2] = 1.0;
+        }
+
+        // Find axis in plane (orthogonal to normal)
+        bool done = false;
+        vtkMath::Cross(normal, axis, px);
+        if (vtkMath::Normalize(px) > 1.0e-3) {
+            done = true;
+        }
+        if (!done) {
+            axis[0] = 0.0;
+            axis[1] = 1.0;
+            axis[2] = 0.0;
+            vtkMath::Cross(normal, axis, px);
+            if (vtkMath::Normalize(px) > 1.0e-3) {
+                done = true;
+            }
+        }
+        if (!done) {
+            axis[0] = 0.0;
+            axis[1] = 0.0;
+            axis[2] = 1.0;
+            vtkMath::Cross(normal, axis, px);
+            vtkMath::Normalize(px);
+        }
+        // Create third orthogonal basis vector
+        vtkMath::Cross(px, normal, py);
+
+        double minSquared = (innerRadius*innerRadius)/(radius*radius);
+        for (int j = 0; j < numPoints; j++) {
+            // Get random sweep angle and radius
+            double angle = getRandomNum(0, 2.0 * vtkMath::DoublePi());
+            // Need sqrt to get uniform distribution
+            double r = sqrt(getRandomNum(minSquared, 1)) * radius;
+            double pt[3];
+            for (int i = 0; i < 3; i++) {
+                pt[i] = center[i] + r * (px[i] * cos(angle) + py[i] * sin(angle));
+            }
+            TRACE("Seed pt: %g %g %g", pt[0], pt[1], pt[2]);
+            pts->InsertNextPoint(pt);
+            cells->InsertNextCell(1);
+            cells->InsertCellPoint(j);
+        }
+
+        seed->SetPoints(pts);
+        seed->SetVerts(cells);
+
+        TRACE("Seed points: %d", seed->GetNumberOfPoints());
+        vtkSmartPointer<vtkDataSet> oldSeed;
+        if (_streamTracer->GetSource() != NULL) {
+            oldSeed = _streamTracer->GetSource();
+        }
+
+        _streamTracer->SetSource(seed);
+        if (oldSeed != NULL) {
+            oldSeed->SetPipelineInformation(NULL);
+        }
+
+        _seedMapper->SetInput(seed);
+    }
+}
+
+/**
  * \brief Use seed points from an n-sided polygon
  *
  * \param[in] center Center point of polygon
  * \param[in] normal Normal vector to orient polygon
+ * \param[in] angle Angle in degrees to rotate about normal
  * \param[in] radius Radius of circumscribing circle
  * \param[in] numSides Number of polygon sides (and points) to generate
@@ -363 +507 @@
 void Streamlines::setSeedToPolygon(double center[3],
                                    double normal[3],
+                                   double angle,
                                    double radius,
                                    int numSides)
@@ -376 +521 @@
         seed->GeneratePolygonOn();
 
+        if (angle != 0.0) {
+            vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
+            trans->RotateWXYZ(angle, normal);
+            vtkSmartPointer<vtkTransformPolyDataFilter> transFilt =
+                vtkSmartPointer<vtkTransformPolyDataFilter>::New();
+            transFilt->SetInputConnection(seed->GetOutputPort());
+            transFilt->SetTransform(trans);
+        }
+
         TRACE("Seed points: %d", numSides);
         vtkSmartPointer<vtkDataSet> oldSeed;
@@ -382 +536 @@
         }
 
-        _streamTracer->SetSourceConnection(seed->GetOutputPort());
+        if (angle != 0.0) {
+            vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
+            trans->Translate(+center[0], +center[1], +center[2]);
+            trans->RotateWXYZ(angle, normal);
+            trans->Translate(-center[0], -center[1], -center[2]);
+            vtkSmartPointer<vtkTransformPolyDataFilter> transFilt =
+                vtkSmartPointer<vtkTransformPolyDataFilter>::New();
+            transFilt->SetInputConnection(seed->GetOutputPort());
+            transFilt->SetTransform(trans);
+            _streamTracer->SetSourceConnection(transFilt->GetOutputPort());
+            _seedMapper->SetInputConnection(transFilt->GetOutputPort());
+        } else {
+            _streamTracer->SetSourceConnection(seed->GetOutputPort());
+            _seedMapper->SetInputConnection(seed->GetOutputPort());
+        }
+
         if (oldSeed != NULL) {
             oldSeed->SetPipelineInformation(NULL);
         }
-
-        _seedMapper->SetInputConnection(seed->GetOutputPort());
+    }
+}
+
+/**
+ * \brief Use seed points from an n-sided polygon
+ *
+ * \param[in] center Center point of polygon
+ * \param[in] normal Normal vector to orient polygon
+ * \param[in] angle Angle in degrees to rotate about normal
+ * \param[in] radius Radius of circumscribing circle
+ * \param[in] numSides Number of polygon sides (and points) to generate
+ * \param[in] numPoints Number of random points to generate
+ */
+void Streamlines::setSeedToFilledPolygon(double center[3],
+                                         double normal[3],
+                                         double angle,
+                                         double radius,
+                                         int numSides,
+                                         int numPoints)
+{
+    if (_streamTracer != NULL) {
+         // Set up seed source object
+        vtkSmartPointer<vtkPolyData> seed = vtkSmartPointer<vtkPolyData>::New();
+        vtkSmartPointer<vtkPoints> pts = vtkSmartPointer<vtkPoints>::New();
+        vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
+
+        // The following code is based on vtkRegularPolygonSource::RequestData
+
+        double px[3];
+        double py[3];
+        double axis[3] = {1., 0., 0.};
+
+        if (vtkMath::Normalize(normal) == 0.0) {
+            normal[0] = 0.0;
+            normal[1] = 0.0;
+            normal[2] = 1.0;
+        }
+
+        // Find axis in plane (orthogonal to normal)
+        bool done = false;
+        vtkMath::Cross(normal, axis, px);
+        if (vtkMath::Normalize(px) > 1.0e-3) {
+            done = true;
+        }
+        if (!done) {
+            axis[0] = 0.0;
+            axis[1] = 1.0;
+            axis[2] = 0.0;
+            vtkMath::Cross(normal, axis, px);
+            if (vtkMath::Normalize(px) > 1.0e-3) {
+                done = true;
+            }
+        }
+        if (!done) {
+            axis[0] = 0.0;
+            axis[1] = 0.0;
+            axis[2] = 1.0;
+            vtkMath::Cross(normal, axis, px);
+            vtkMath::Normalize(px);
+        }
+        // Create third orthogonal basis vector
+        vtkMath::Cross(px, normal, py);
+
+        double verts[numSides][3];
+        double sliceTheta = 2.0 * vtkMath::DoublePi() / (double)numSides;
+        angle = vtkMath::RadiansFromDegrees(angle);
+        for (int j = 0; j < numSides; j++) {
+            for (int i = 0; i < 3; i++) {
+                double theta = sliceTheta * (double)j - angle;
+                verts[j][i] = center[i] + radius * (px[i] * cos(theta) + 
+                                                    py[i] * sin(theta));
+            }
+            TRACE("Vert %d: %g %g %g", j, verts[j][0], verts[j][1], verts[j][2]);
+        }
+
+        // Note: this gives a uniform distribution because the polygon is regular and
+        // the triangular sections have equal area
+        if (numSides == 3) {
+            for (int j = 0; j < numPoints; j++) {
+                double pt[3];
+                getRandomPointInTriangle(pt, verts[0], verts[1], verts[2]);
+                TRACE("Seed pt: %g %g %g", pt[0], pt[1], pt[2]);
+                pts->InsertNextPoint(pt);
+                cells->InsertNextCell(1);
+                cells->InsertCellPoint(j);
+            }
+        } else {
+            for (int j = 0; j < numPoints; j++) {
+                // Get random triangle section
+                int tri = rand() % numSides;
+                double pt[3];
+                getRandomPointInTriangle(pt, center, verts[tri], verts[(tri+1) % numSides]);
+                TRACE("Seed pt: %g %g %g", pt[0], pt[1], pt[2]);
+                pts->InsertNextPoint(pt);
+                cells->InsertNextCell(1);
+                cells->InsertCellPoint(j);
+            }
+        }
+
+        seed->SetPoints(pts);
+        seed->SetVerts(cells);
+
+        TRACE("Seed points: %d", seed->GetNumberOfPoints());
+        vtkSmartPointer<vtkDataSet> oldSeed;
+        if (_streamTracer->GetSource() != NULL) {
+            oldSeed = _streamTracer->GetSource();
+        }
+
+        _streamTracer->SetSource(seed);
+        if (oldSeed != NULL) {
+            oldSeed->SetPipelineInformation(NULL);
+        }
+
+        _seedMapper->SetInput(seed);
     }
 }

trunk/packages/vizservers/vtkvis/RpStreamlines.h

@@ -64 +64 @@
     void setSeedToRake(double start[3], double end[3], int numPoints);
 
+    void setSeedToDisk(double center[3], double normal[3],
+                       double radius, double innerRadius, int numPoints);
+
     void setSeedToPolygon(double center[3], double normal[3],
-                          double radius, int numSides);
+                          double angle, double radius,
+                          int numSides);
+
+    void setSeedToFilledPolygon(double center[3], double normal[3],
+                                double angle, double radius,
+                                int numSides, int numPoints);
 
     void setMaxPropagation(double length);
@@ -87 +95 @@
     virtual void update();
 
+    static double getRandomNum(double min, double max);
     static void getRandomPoint(double pt[3], const double bounds[6]);
-    static void getRandomCellPt(vtkDataSet *ds, double pt[3]);
+    static void getRandomPointInTriangle(double pt[3],
+                                         const double v1[3],
+                                         const double v2[3],
+                                         const double v3[3]);
+    static void getRandomPointOnLineSegment(double pt[3],
+                                            const double endpt[3],
+                                            const double endpt2[3]);
+    static void getRandomCellPt(double pt[3], vtkDataSet *ds);
 
     LineType _lineType;

trunk/packages/vizservers/vtkvis/RpVtkRenderer.cpp

@@ -5442 +5442 @@
  * \brief Set the streamlines seed to points along a line
  */
-void Renderer::setStreamlinesSeedToRake(const DataSetId& id, double start[3], double end[3], int numPoints)
+void Renderer::setStreamlinesSeedToRake(const DataSetId& id,
+                                        double start[3], double end[3],
+                                        int numPoints)
 {
     StreamlinesHashmap::iterator itr;
@@ -5467 +5469 @@
 
 /**
+ * \brief Set the streamlines seed to points inside a disk, with optional
+ * hole
+ */
+void Renderer::setStreamlinesSeedToDisk(const DataSetId& id,
+                                        double center[3], double normal[3],
+                                        double radius, double innerRadius,
+                                        int numPoints)
+{
+    StreamlinesHashmap::iterator itr;
+
+    bool doAll = false;
+
+    if (id.compare("all") == 0) {
+        itr = _streamlines.begin();
+        doAll = true;
+    } else {
+        itr = _streamlines.find(id);
+    }
+    if (itr == _streamlines.end()) {
+        ERROR("Streamlines not found: %s", id.c_str());
+        return;
+    }
+
+    do {
+        itr->second->setSeedToDisk(center, normal, radius, innerRadius, numPoints);
+    } while (doAll && ++itr != _streamlines.end());
+
+    _needsRedraw = true;
+}
+
+/**
  * \brief Set the streamlines seed to vertices of an n-sided polygon
  */
 void Renderer::setStreamlinesSeedToPolygon(const DataSetId& id,
                                            double center[3], double normal[3],
-                                           double radius, int numSides)
+                                           double angle, double radius,
+                                           int numSides)
 {
     StreamlinesHashmap::iterator itr;
@@ -5489 +5523 @@
 
     do {
-        itr->second->setSeedToPolygon(center, normal, radius, numSides);
+        itr->second->setSeedToPolygon(center, normal, angle, radius, numSides);
+    } while (doAll && ++itr != _streamlines.end());
+
+    _needsRedraw = true;
+}
+
+/**
+ * \brief Set the streamlines seed to vertices of an n-sided polygon
+ */
+void Renderer::setStreamlinesSeedToFilledPolygon(const DataSetId& id,
+                                                 double center[3],
+                                                 double normal[3],
+                                                 double angle, double radius,
+                                                 int numSides, int numPoints)
+{
+    StreamlinesHashmap::iterator itr;
+
+    bool doAll = false;
+
+    if (id.compare("all") == 0) {
+        itr = _streamlines.begin();
+        doAll = true;
+    } else {
+        itr = _streamlines.find(id);
+    }
+    if (itr == _streamlines.end()) {
+        ERROR("Streamlines not found: %s", id.c_str());
+        return;
+    }
+
+    do {
+        itr->second->setSeedToFilledPolygon(center, normal, angle,
+                                            radius, numSides, numPoints);
     } while (doAll && ++itr != _streamlines.end());
 
@@ -6388 +6454 @@
 
 /**
+ * \brief Set the VTK camera parameters based on a 4x4 view matrix
+ */
+void Renderer::setCameraFromMatrix(vtkCamera *camera, vtkMatrix4x4& mat)
+{
+    double d = camera->GetDistance();
+    double vu[3];
+    vu[0] = mat[1][0];
+    vu[1] = mat[1][1];
+    vu[2] = mat[1][2];
+    double trans[3];
+    trans[0] = mat[0][3];
+    trans[1] = mat[1][3];
+    trans[2] = mat[2][3];
+    mat[0][3] = 0;
+    mat[1][3] = 0;
+    mat[2][3] = 0;
+    double vpn[3] = {mat[2][0], mat[2][1], mat[2][2]};
+    double pos[3];
+    // With translation removed, we have an orthogonal matrix, 
+    // so the inverse is the transpose
+    mat.Transpose();
+    mat.MultiplyPoint(trans, pos);
+    pos[0] = -pos[0];
+    pos[1] = -pos[1];
+    pos[2] = -pos[2];
+    double fp[3];
+    fp[0] = pos[0] - vpn[0] * d;
+    fp[1] = pos[1] - vpn[1] * d;
+    fp[2] = pos[2] - vpn[2] * d;
+    camera->SetPosition(pos);
+    camera->SetFocalPoint(fp);
+    camera->SetViewUp(vu);
+}
+
+inline void quaternionToMatrix4x4(double quat[4], vtkMatrix4x4& mat)
+{
+    double ww = quat[0]*quat[0];
+    double wx = quat[0]*quat[1];
+    double wy = quat[0]*quat[2];
+    double wz = quat[0]*quat[3];
+
+    double xx = quat[1]*quat[1];
+    double yy = quat[2]*quat[2];
+    double zz = quat[3]*quat[3];
+
+    double xy = quat[1]*quat[2];
+    double xz = quat[1]*quat[3];
+    double yz = quat[2]*quat[3];
+
+    double rr = xx + yy + zz;
+    // normalization factor, just in case quaternion was not normalized
+    double f = double(1)/double(sqrt(ww + rr));
+    double s = (ww - rr)*f;
+    f *= 2;
+
+    mat[0][0] = xx*f + s;
+    mat[1][0] = (xy + wz)*f;
+    mat[2][0] = (xz - wy)*f;
+
+    mat[0][1] = (xy - wz)*f;
+    mat[1][1] = yy*f + s;
+    mat[2][1] = (yz + wx)*f;
+
+    mat[0][2] = (xz + wy)*f;
+    mat[1][2] = (yz - wx)*f;
+    mat[2][2] = zz*f + s;
+}
+
+inline void quaternionToTransposeMatrix4x4(double quat[4], vtkMatrix4x4& mat)
+{
+    double ww = quat[0]*quat[0];
+    double wx = quat[0]*quat[1];
+    double wy = quat[0]*quat[2];
+    double wz = quat[0]*quat[3];
+
+    double xx = quat[1]*quat[1];
+    double yy = quat[2]*quat[2];
+    double zz = quat[3]*quat[3];
+
+    double xy = quat[1]*quat[2];
+    double xz = quat[1]*quat[3];
+    double yz = quat[2]*quat[3];
+
+    double rr = xx + yy + zz;
+    // normalization factor, just in case quaternion was not normalized
+    double f = double(1)/double(sqrt(ww + rr));
+    double s = (ww - rr)*f;
+    f *= 2;
+
+    mat[0][0] = xx*f + s;
+    mat[0][1] = (xy + wz)*f;
+    mat[0][2] = (xz - wy)*f;
+
+    mat[1][0] = (xy - wz)*f;
+    mat[1][1] = yy*f + s;
+    mat[1][2] = (yz + wx)*f;
+
+    mat[2][0] = (xz + wy)*f;
+    mat[2][1] = (yz - wx)*f;
+    mat[2][2] = zz*f + s;
+}
+
+/**
  * \brief Set the orientation of the camera from a quaternion
  *
  * \param[in] quat A quaternion with scalar part first: w,x,y,z
- */
-void Renderer::setCameraOrientation(double quat[4])
+ * \param[in] absolute Is rotation absolute or relative?
+ */
+void Renderer::setCameraOrientation(double quat[4], bool absolute)
 {
     if (_cameraMode == IMAGE)
@@ -6398 +6568 @@
     vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
     vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
-    double mat3[3][3];
-    vtkMath::QuaternionToMatrix3x3(quat, mat3);
     vtkSmartPointer<vtkMatrix4x4> mat4 = vtkSmartPointer<vtkMatrix4x4>::New();
-    for (int r = 0; r < 3; r++) {
-        memcpy((*mat4)[r], mat3[r], sizeof(double)*3);
-    }
-    TRACE("Arcball camera matrix: %g %g %g %g %g %g %g %g %g",
-          (*mat4)[0][0], (*mat4)[0][1], (*mat4)[0][2],
-          (*mat4)[1][0], (*mat4)[1][1], (*mat4)[1][2],
-          (*mat4)[2][0], (*mat4)[2][1], (*mat4)[2][2]);
-    camera->SetPosition(0, 0, 1);
-    camera->SetFocalPoint(0, 0, 0);
-    camera->SetViewUp(0, 1, 0);
-    setViewAngle(_windowHeight);
-    double bounds[6];
-    collectBounds(bounds, false);
-    _renderer->ResetCamera(bounds);
-    camera->GetFocalPoint(_cameraFocalPoint);
-    trans->Translate(+_cameraFocalPoint[0], +_cameraFocalPoint[1], +_cameraFocalPoint[2]);
-    trans->Concatenate(mat4);
-    trans->Translate(-_cameraFocalPoint[0], -_cameraFocalPoint[1], -_cameraFocalPoint[2]);
-    camera->ApplyTransform(trans);
+
+    if (absolute) {
+        quaternionToMatrix4x4(quat, *mat4);
+#ifdef DEBUG
+        TRACE("Arcball camera matrix:\n %g %g %g\n %g %g %g\n %g %g %g",
+              (*mat4)[0][0], (*mat4)[0][1], (*mat4)[0][2],
+              (*mat4)[1][0], (*mat4)[1][1], (*mat4)[1][2],
+              (*mat4)[2][0], (*mat4)[2][1], (*mat4)[2][2]);
+#endif
+        camera->SetPosition(0, 0, 1);
+        camera->SetFocalPoint(0, 0, 0);
+        camera->SetViewUp(0, 1, 0);
+        double bounds[6];
+        collectBounds(bounds, false);
+        _renderer->ResetCamera(bounds);
+        camera->GetFocalPoint(_cameraFocalPoint);
+        trans->Translate(+_cameraFocalPoint[0], +_cameraFocalPoint[1], +_cameraFocalPoint[2]);
+        trans->Concatenate(mat4);
+        trans->Translate(-_cameraFocalPoint[0], -_cameraFocalPoint[1], -_cameraFocalPoint[2]);
+        camera->ApplyTransform(trans);
+    } else {
+        quaternionToTransposeMatrix4x4(quat, *mat4);
+#ifdef DEBUG
+        vtkSmartPointer<vtkMatrix4x4> mat2 = vtkSmartPointer<vtkMatrix4x4>::New();
+        mat2->DeepCopy(camera->GetViewTransformMatrix());
+        TRACE("camera matrix:\n %g %g %g %g\n %g %g %g %g\n %g %g %g %g\n %g %g %g %g",
+              (*mat2)[0][0], (*mat2)[0][1], (*mat2)[0][2], (*mat2)[0][3],
+              (*mat2)[1][0], (*mat2)[1][1], (*mat2)[1][2], (*mat2)[1][3],
+              (*mat2)[2][0], (*mat2)[2][1], (*mat2)[2][2], (*mat2)[2][3],
+              (*mat2)[3][0], (*mat2)[3][1], (*mat2)[3][2], (*mat2)[3][3]);
+        printCameraInfo(camera);
+#endif
+        trans->Translate(0, 0, -camera->GetDistance());
+        trans->Concatenate(mat4);
+        trans->Translate(0, 0, camera->GetDistance());
+        trans->Concatenate(camera->GetViewTransformMatrix());
+        setCameraFromMatrix(camera, *trans->GetMatrix());
+    }
     storeCameraOrientation();
-    if (_cameraZoomRatio != 1.0) {
-        double z = _cameraZoomRatio;
-        _cameraZoomRatio = 1.0;
-        zoomCamera(z, true);
-    }
-    if (_cameraPan[0] != 0.0 || _cameraPan[1] != 0.0) {
-        double panx = _cameraPan[0];
-        double pany = -_cameraPan[1];
-        _cameraPan[0] = 0;
-        _cameraPan[1] = 0;
-        panCamera(panx, pany, true);
+    if (absolute) {
+        if (_cameraZoomRatio != 1.0) {
+            double z = _cameraZoomRatio;
+            _cameraZoomRatio = 1.0;
+            zoomCamera(z, true);
+        }
+        if (_cameraPan[0] != 0.0 || _cameraPan[1] != 0.0) {
+            double panx = _cameraPan[0];
+            double pany = -_cameraPan[1];
+            _cameraPan[0] = 0;
+            _cameraPan[1] = 0;
+            panCamera(panx, pany, true);
+        }
     }
     _renderer->ResetCameraClippingRange();
@@ -6521 +6710 @@
         _renderer->ResetCamera(bounds);
         _renderer->ResetCameraClippingRange();
-        computeScreenWorldCoords();
+        //computeScreenWorldCoords();
     }
 
@@ -6550 +6739 @@
     _renderer->ResetCameraClippingRange();
     storeCameraOrientation();
-    computeScreenWorldCoords();
+    //computeScreenWorldCoords();
     _needsRedraw = true;
 }
@@ -7193 +7382 @@
         resetAxes();
         _renderer->ResetCamera(bounds);
-        computeScreenWorldCoords();
+        //computeScreenWorldCoords();
         break;
     case PERSPECTIVE:
@@ -7199 +7388 @@
         resetAxes();
         _renderer->ResetCamera(bounds);
-        computeScreenWorldCoords();
+        //computeScreenWorldCoords();
         break;
     default:

trunk/packages/vizservers/vtkvis/RpVtkRenderer.h

@@ -150 +150 @@
     void rotateCamera(double yaw, double pitch, double roll);
 
-    void setCameraOrientation(double quat[4]);
+    void setCameraOrientation(double quat[4], bool absolute = true);
+
+    void panCamera(double x, double y, bool absolute = true);
+
+    void zoomCamera(double z, bool absolute = true);
 
     void setCameraOrientationAndPosition(double position[3],
@@ -159 +163 @@
                                          double focalPoint[3],
                                          double viewUp[3]);
-
-    void panCamera(double x, double y, bool absolute = true);
-
-    void zoomCamera(double z, bool absolute = true);
 
     // Rendering an image
@@ -559 +559 @@
     void setStreamlinesSeedToRandomPoints(const DataSetId& id, int numPoints);
 
-    void setStreamlinesSeedToRake(const DataSetId& id, double start[3], double end[3], int numPoints);
-
-    void setStreamlinesSeedToPolygon(const DataSetId& id, double center[3], double normal[3],
-                                           double radius, int numSides);
+    void setStreamlinesSeedToRake(const DataSetId& id,
+                                  double start[3], double end[3],
+                                  int numPoints);
+
+    void setStreamlinesSeedToDisk(const DataSetId& id,
+                                  double center[3], double normal[3],
+                                  double radius, double innerRadius,
+                                  int numPoints);
+
+    void setStreamlinesSeedToPolygon(const DataSetId& id,
+                                     double center[3], double normal[3],
+                                     double angle, double radius,
+                                     int numSides);
+
+    void setStreamlinesSeedToFilledPolygon(const DataSetId& id,
+                                           double center[3], double normal[3],
+                                           double angle, double radius,
+                                           int numSides, int numPoints);
 
     void setStreamlinesLength(const DataSetId& id, double length);
@@ -627 +641 @@
 
     void updateRanges(bool useCumulative);
+
+    void setCameraFromMatrix(vtkCamera *camera, vtkMatrix4x4 &mat);
 
     void computeDisplayToWorld(double x, double y, double z, double worldPt[4]);

trunk/packages/vizservers/vtkvis/RpVtkRendererCmd.cpp

@@ -3720 +3720 @@
 
 static int
+StreamlinesSeedDiskOp(ClientData clientData, Tcl_Interp *interp, int objc, 
+                      Tcl_Obj *const *objv)
+{
+    double center[3], normal[3], radius, innerRadius;
+    int numPoints;
+    for (int i = 0; i < 3; i++) {
+        if (Tcl_GetDoubleFromObj(interp, objv[3+i], &center[i]) != TCL_OK) {
+            return TCL_ERROR;
+        }
+        if (Tcl_GetDoubleFromObj(interp, objv[6+i], &normal[i]) != TCL_OK) {
+            return TCL_ERROR;
+        }
+    }
+    if (Tcl_GetDoubleFromObj(interp, objv[9], &radius) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (Tcl_GetDoubleFromObj(interp, objv[10], &innerRadius) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (Tcl_GetIntFromObj(interp, objv[11], &numPoints) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (objc == 13) {
+        const char *name = Tcl_GetString(objv[12]);
+        g_renderer->setStreamlinesSeedToDisk(name, center, normal, radius, innerRadius, numPoints);
+    } else {
+        g_renderer->setStreamlinesSeedToDisk("all", center, normal, radius, innerRadius, numPoints);
+    }
+    return TCL_OK;
+}
+
+static int
 StreamlinesSeedPolygonOp(ClientData clientData, Tcl_Interp *interp, int objc, 
                          Tcl_Obj *const *objv)
 {
-    double center[3], normal[3], radius;
+    double center[3], normal[3], angle, radius;
     int numSides;
     for (int i = 0; i < 3; i++) {
@@ -3733 +3765 @@
         }
     }
-    if (Tcl_GetDoubleFromObj(interp, objv[9], &radius) != TCL_OK) {
-        return TCL_ERROR;
-    }
-    if (Tcl_GetIntFromObj(interp, objv[10], &numSides) != TCL_OK) {
-        return TCL_ERROR;
-    }
-    if (objc == 12) {
-        const char *name = Tcl_GetString(objv[11]);
-        g_renderer->setStreamlinesSeedToPolygon(name, center, normal, radius, numSides);
-    } else {
-        g_renderer->setStreamlinesSeedToPolygon("all", center, normal, radius, numSides);
+    if (Tcl_GetDoubleFromObj(interp, objv[9], &angle) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (Tcl_GetDoubleFromObj(interp, objv[10], &radius) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (Tcl_GetIntFromObj(interp, objv[11], &numSides) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (numSides < 3) {
+        Tcl_AppendResult(interp, "numSides must be 3 or greater", (char*)NULL);
+        return TCL_ERROR;
+    }
+    if (objc == 13) {
+        const char *name = Tcl_GetString(objv[12]);
+        g_renderer->setStreamlinesSeedToPolygon(name, center, normal, angle, radius, numSides);
+    } else {
+        g_renderer->setStreamlinesSeedToPolygon("all", center, normal, angle, radius, numSides);
+    }
+    return TCL_OK;
+}
+
+static int
+StreamlinesSeedFilledPolygonOp(ClientData clientData, Tcl_Interp *interp, int objc, 
+                               Tcl_Obj *const *objv)
+{
+    double center[3], normal[3], angle, radius;
+    int numSides, numPoints;
+    for (int i = 0; i < 3; i++) {
+        if (Tcl_GetDoubleFromObj(interp, objv[3+i], &center[i]) != TCL_OK) {
+            return TCL_ERROR;
+        }
+        if (Tcl_GetDoubleFromObj(interp, objv[6+i], &normal[i]) != TCL_OK) {
+            return TCL_ERROR;
+        }
+    }
+    if (Tcl_GetDoubleFromObj(interp, objv[9], &angle) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (Tcl_GetDoubleFromObj(interp, objv[10], &radius) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (Tcl_GetIntFromObj(interp, objv[11], &numSides) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (numSides < 3) {
+        Tcl_AppendResult(interp, "numSides must be 3 or greater", (char*)NULL);
+        return TCL_ERROR;
+    }
+    if (Tcl_GetIntFromObj(interp, objv[12], &numPoints) != TCL_OK) {
+        return TCL_ERROR;
+    }
+    if (objc == 14) {
+        const char *name = Tcl_GetString(objv[13]);
+        g_renderer->setStreamlinesSeedToFilledPolygon(name, center, normal, angle,
+                                                      radius, numSides, numPoints);
+    } else {
+        g_renderer->setStreamlinesSeedToFilledPolygon("all", center, normal, angle,
+                                                      radius, numSides, numPoints);
     }
     return TCL_OK;
@@ -3810 +3890 @@
 static Rappture::CmdSpec streamlinesSeedOps[] = {
     {"color",   1, StreamlinesSeedColorOp, 6, 7, "r g b ?dataSetName?"},
-    {"polygon", 1, StreamlinesSeedPolygonOp, 11, 12, "centerX centerY centerZ normalX normalY normalZ radius numSides ?dataSetName?"},
+    {"disk",    1, StreamlinesSeedDiskOp, 12, 13, "centerX centerY centerZ normalX normalY normalZ radius innerRadius numPoints ?dataSetName?"},
+    {"fpoly",   1, StreamlinesSeedFilledPolygonOp, 13, 14, "centerX centerY centerZ normalX normalY normalZ angle radius numSides numPoints ?dataSetName?"},
+    {"polygon", 1, StreamlinesSeedPolygonOp, 12, 13, "centerX centerY centerZ normalX normalY normalZ angle radius numSides ?dataSetName?"},
     {"rake",    3, StreamlinesSeedRakeOp, 10, 11, "startX startY startZ endX endY endZ numPoints ?dataSetName?"},
     {"random",  3, StreamlinesSeedRandomOp, 4, 5, "numPoints ?dataSetName?"},
@@ -3884 +3966 @@
     {"ribbons",   1, StreamlinesRibbonsOp, 4, 5, "width angle ?dataSetName?"},
     {"scale",     2, StreamlinesScaleOp, 5, 6, "sx sy sz ?dataSetName?"},
-    {"seed",      2, StreamlinesSeedOp, 4, 11, "op params... ?dataSetName?"},
+    {"seed",      2, StreamlinesSeedOp, 4, 14, "op params... ?dataSetName?"},
     {"tubes",     1, StreamlinesTubesOp, 4, 5, "numSides radius ?dataSetName?"},
     {"visible",   1, StreamlinesVisibleOp, 3, 4, "bool ?dataSetName?"}

trunk/packages/vizservers/vtkvis/protocol.txt

@@ -55 +55 @@
        Option all resets orientation/rotation as well as pan/zoom/clip range
 camera rotate <yaw> <pitch> <roll>
-       Specify relative rotation in Euler angles (FIXME)
+       Specify relative rotation in Euler angles
 camera set <posX> <posY> <posZ> <focalPtX> <focalPtY> <focalPtZ> <viewUpX> <viewUpY> <viewUpZ>
        Set camera orientation using camera position, focal point and view up 
@@ -250 +250 @@
 streamlines scale <sx> <sy> <sz> <?dataSetName?>
 streamlines seed color <r> <g> <b> <?datasetName?>
-streamlines seed polygon <centerX> <centerY> <centerZ> <normalX> <normalY> <normalZ> <radius> <numSides> <?dataSetName?>
+streamlines seed disk <centerX> <centerY> <centerZ> <normalX> <normalY> <normalZ> <radius> <innerRadius> <numPoints> <?dataSetName?>
+            Create a disk seed area with optional hole, filled with randomly placed points
+streamlines seed fpoly <centerX> <centerY> <centerZ> <normalX> <normalY> <normalZ> <angle> <radius> <numSides> <numPoints> <?dataSetName?>
+            Create a regular n-sided polygonal seed area filled with randomly placed points
+streamlines seed polygon <centerX> <centerY> <centerZ> <normalX> <normalY> <normalZ> <angle> <radius> <numSides> <?dataSetName?>
+            Create seed points from vertices of a regular n-sided polygon
 streamlines seed rake <startX> <startY> <startZ> <endX> <endY> <endZ> <numPoints> <?datasetName?>
 streamlines seed random <numPoints> <?datasetName?>