/Main/src/Microsoft.Research.Visualization3D/Providers/Isosurfaces/Metaball.cs

  1using System;
  2using System.Collections.Generic;
  3using Microsoft.Research.Visualization3D.VertexStructures;
  4using SlimDX.Direct3D9;
  5using Microsoft.Research.Visualization3D.CameraUtilities;
  6using SlimDX;
  7using Microsoft.Research.Visualization3D.Auxilaries;
  8using System.Windows.Threading;
  9using System.Diagnostics;
 10using Microsoft.Research.Visualization3D.MainLoops;
 11using System.Threading;
 12
 13namespace Microsoft.Research.Visualization3D.Isosurfaces
 14{
 15    public partial class Metaball : DrawableComponent
 16    {
 17        List<MetaballVertex> vertexList;
 18        int primitiveCount;
 19
 20        int iDataSetSizeX;
 21        int iDataSetSizeY;
 22        int iDataSetSizeZ;
 23
 24        VertexBuffer vb;
 25
 26        public Metaball(DX3DHost host, Visualization3DDataSource dataSource) :
 27            base(host, dataSource)
 28        {
 29            iDataSetSizeX = dataSource.DisplayData.GetLength(0);
 30            iDataSetSizeY = dataSource.DisplayData.GetLength(1);
 31            iDataSetSizeZ = dataSource.DisplayData.GetLength(2);
 32            this.WpfDispatcher = Dispatcher.CurrentDispatcher;
 33        }
 34
 35        public override void Initialize()
 36        {
 37            
 38            if (effect == null)
 39            {
 40                effect = Effect.FromStream(device, this.GetType().Assembly.GetManifestResourceStream("Microsoft.Research.Visualization3D.Shaders.PerPixelLightning.fx"), ShaderFlags.None);
 41            }
 42            SetCamera();
 43            vertexList = new List<MetaballVertex>();
 44            Color3 mbColor = RgbPalette.GetColor(fTargetValue, dataSource.Maximum, dataSource.Minimum, dataSource.MissingValue);
 45            effect.SetValue("mbColor", new Vector4(mbColor.Red, mbColor.Green, mbColor.Blue, 1.0f));
 46
 47            MarchingCubes();
 48            Completed();
 49
 50            base.Initialize();
 51        }
 52
 53        //Calling Cube building function for each voxel
 54        public void MarchingCubes()
 55        {
 56            vertexList = new List<MetaballVertex>();
 57
 58            for (int iX = 0; iX < iDataSetSizeX - 1; iX++)
 59                for (int iY = 0; iY < iDataSetSizeY - 1; iY++)
 60                    for (int iZ = 0; iZ < iDataSetSizeZ - 1; iZ++)
 61                    {
 62                        vMarchCube(iX, iY, iZ, 1.0f);
 63                    }
 64        }
 65
 66        float fGetOffset(float fValue1, float fValue2, float fValueDesired)
 67        {
 68            double fDelta = fValue2 - fValue1;
 69
 70            if (fDelta == 0.0)
 71            {
 72                return 0.5f;
 73            }
 74            return (fValueDesired - fValue1) / (float)fDelta;
 75        }
 76
 77        float fSample(float fX, float fY, float fZ)
 78        {
 79            return MathHelper.GetValue(new Vector3(fX, fY, fZ), dataSource.DisplayData);
 80        }
 81
 82        //vGetColor fins color via point and normal
 83        void vGetColor(ref Vector3 rfColor, Vector3 rfPosition, Vector3 rfNormal)
 84        {
 85            float fX = rfNormal.X;
 86            float fY = rfNormal.Y;
 87            float fZ = rfNormal.Z;
 88            rfColor.X = (float)((fX > 0.0 ? fX : 0.0) + (fY < 0.0 ? -0.5 * fY : 0.0) + (fZ < 0.0 ? -0.5 * fZ : 0.0));
 89            rfColor.Y = (float)((fY > 0.0 ? fY : 0.0) + (fZ < 0.0 ? -0.5 * fZ : 0.0) + (fX < 0.0 ? -0.5 * fX : 0.0));
 90            rfColor.Z = (float)((fZ > 0.0 ? fZ : 0.0) + (fX < 0.0 ? -0.5 * fX : 0.0) + (fY < 0.0 ? -0.5 * fY : 0.0));
 91        }
 92
 93        //vGetNormal() Calculates normal via gradient to point
 94        void vGetNormal(ref Vector3 rfNormal, float fX, float fY, float fZ)
 95        {
 96            rfNormal.X = MathHelper.GetValue(new Vector3(fX - 0.01f, fY, fZ), dataSource.DisplayData) - MathHelper.GetValue(new Vector3(fX + 0.01f, fY, fZ), dataSource.DisplayData);
 97            rfNormal.Y = MathHelper.GetValue(new Vector3(fX, fY - 0.01f, fZ), dataSource.DisplayData) - MathHelper.GetValue(new Vector3(fX, fY + 0.01f, fZ), dataSource.DisplayData);
 98            rfNormal.Z = MathHelper.GetValue(new Vector3(fX, fY, fZ - 0.01f), dataSource.DisplayData) - MathHelper.GetValue(new Vector3(fX, fY, fZ + 0.01f), dataSource.DisplayData);
 99            if (rfNormal.Length() > 0)
100                rfNormal.Normalize();
101        }
102
103        //Find part of Surface for current voxel
104        void vMarchCube(float fX, float fY, float fZ, float fScale)
105        {
106
107
108            int iCorner, iVertex, iVertexTest, iEdge, iTriangle, iFlagIndex, iEdgeFlags;
109            float fOffset;
110            Vector3 sColor = Vector3.Zero;
111            float[] afCubeValue = new float[8];
112            Vector3[] asEdgeVertex = new Vector3[12];
113            Vector3[] asEdgeNorm = new Vector3[12];
114
115            //Find value in voxel's knots
116            for (iVertex = 0; iVertex < 8; iVertex++)
117            {
118                afCubeValue[iVertex] = (float)dataSource.DisplayData[(int)(fX + a2fVertexOffset[iVertex, 0]), (int)(fY + a2fVertexOffset[iVertex, 1]), (int)(fZ + a2fVertexOffset[iVertex, 2])];
119                if (afCubeValue[iVertex] == dataSource.MissingValue) return;
120            }
121
122            //Checking intersections via table
123            iFlagIndex = 0;
124            for (iVertexTest = 0; iVertexTest < 8; iVertexTest++)
125            {
126                if (afCubeValue[iVertexTest] <= fTargetValue)
127                    iFlagIndex |= 1 << iVertexTest;
128            }
129
130            //Finally, get our surface configuration
131            iEdgeFlags = aiCubeEdgeFlags[iFlagIndex];
132
133            //if our surface doesn't intersect voxel
134            if (iEdgeFlags == 0)
135            {
136                return;
137            }
138
139            //Building vertices
140            for (iEdge = 0; iEdge < 12; iEdge++)
141            {
142                if ((iEdgeFlags & (1 << iEdge)) != 0)
143                {
144                    fOffset = fGetOffset(afCubeValue[a2iEdgeConnection[iEdge, 0]],
145                               afCubeValue[a2iEdgeConnection[iEdge, 1]], fTargetValue);
146
147                    asEdgeVertex[iEdge].X = (float)(fX + (a2fVertexOffset[a2iEdgeConnection[iEdge, 0], 0] + fOffset * a2fEdgeDirection[iEdge, 0]) * fScale);
148                    asEdgeVertex[iEdge].Y = (float)(fY + (a2fVertexOffset[a2iEdgeConnection[iEdge, 0], 1] + fOffset * a2fEdgeDirection[iEdge, 1]) * fScale);
149                    asEdgeVertex[iEdge].Z = (float)(fZ + (a2fVertexOffset[a2iEdgeConnection[iEdge, 0], 2] + fOffset * a2fEdgeDirection[iEdge, 2]) * fScale);
150
151                    vGetNormal(ref asEdgeNorm[iEdge], asEdgeVertex[iEdge].X, asEdgeVertex[iEdge].Y, asEdgeVertex[iEdge].Z);
152                }
153            }
154
155
156            //Building triangles
157            for (iTriangle = 0; iTriangle < 5; iTriangle++)
158            {
159                if (a2iTriangleConnectionTable[iFlagIndex, 3 * iTriangle] < 0)
160                    break;
161
162                for (iCorner = 0; iCorner < 3; iCorner++)
163                {
164                    iVertex = a2iTriangleConnectionTable[iFlagIndex, 3 * iTriangle + iCorner];
165
166                    vertexList.Add(new MetaballVertex(
167                        asEdgeVertex[iVertex],
168                        asEdgeNorm[iVertex]
169                        ));
170
171                }
172
173            }
174
175
176        }
177
178        protected override void SetCamera()
179        {
180            camera.CameraScale = 12.0f;
181            
182            float cameraScale = 2.0f;
183            camera.Location = cameraScale * new Vector3(dataSource.DisplayData.GetLength(0), dataSource.DisplayData.GetLength(1), dataSource.DisplayData.GetLength(2));
184            camera.Target = new Vector3(dataSource.DisplayData.GetLength(0) / 2f, dataSource.DisplayData.GetLength(1) / 2f, dataSource.DisplayData.GetLength(2) / 2f);
185            camera.Up = new Vector3(0, 1, 0);
186
187            effect.SetValue("world", Matrix.Identity);
188            effect.SetValue("view", camera.ViewMatrix);
189            effect.SetValue("projection", camera.ProjectionMatrix);
190
191            effect.SetValue("cameraPosition", camera.Location);
192            effect.SetValue("lightPosition", camera.Location);
193            effect.SetValue("ambientLightColor", new Vector4(1.0f, 1.0f, 1.0f, 1));
194            effect.SetValue("diffuseLightColor", new Vector4(1.0f, 1.0f, 1.0f, 1));
195            effect.SetValue("specularLightColor", new Vector4(1.0f, 1.0f, 1.0f, 1));
196
197            effect.SetValue("specularPower", 1.0f);
198            effect.SetValue("specularIntensity", 1.0f);
199        }
200
201        public override void Draw(TimeEntity timeEntity)
202        {
203            if (primitiveCount > 0)
204            {
205                effect.Technique = new EffectHandle("PerPixelDiffuseAndPhongMetaball");
206                device.SetStreamSource(0, vb, 0, MetaballVertex.SizeInBytes);
207                device.VertexFormat = MetaballVertex.Format;
208                int passes = effect.Begin();
209                for (int i = 0; i < passes; i++)
210                {
211                    effect.BeginPass(i);
212                    device.DrawPrimitives(PrimitiveType.TriangleList, 0, primitiveCount);
213                    effect.EndPass();
214                }
215                effect.End();
216            }
217        }
218
219        public override void Update(TimeEntity timeEntity)
220        {
221            effect.SetValue("world", Matrix.Identity);
222            effect.SetValue("view", camera.ViewMatrix);
223            effect.SetValue("projection", camera.ProjectionMatrix);
224
225            effect.SetValue("cameraPosition", camera.Location);
226            effect.SetValue("lightPosition", camera.Location);
227        }
228    }
229}