feat(metrics): add BDR, RaR, VR, SVF and AR urban form indices

CHANGELOG.md

@@ -6,6 +6,13 @@ and this project adheres to [Semantic Versioning](https://semver.org/).
 
 ## [Unreleased]
 
+### Added
+- Urban form indices in the morphology metrics: built density ratio (BDR), road area ratio (RaR), vertical area ratio (VR), mean ground-level sky view factor (SVF) and aspect ratio (AR), with the underlying built volume, road area and facade area. SVF is ray-cast from a grid of roughly 400 open-ground points, 144 cosine-weighted rays each, against the building masses and the terrain when one is wired; trees do not occlude. AR is the mean building height over the `StreetWidth` input. Road area is the site minus every parcel, parks included.
+- The existing plan area density and plan-area-weighted mean height are now labelled BCR and BHt in the `MorphologyMetrics` text; their manifest keys are unchanged.
+
+### Changed
+- Manifest schema version is now `1.2.0` for the added morphology keys. The schema version is part of the case identifier, so identifiers differ from those produced by 0.1.0.5.
+
 ## [0.1.0.5] - 2026-08-27
 
 ### Added

README.md

@@ -75,7 +75,7 @@ Right-click the Massing Generator and use **Street Network** to select `Irregula
 
 **Outputs**: Footprints, Masses, Heights, Streets, FloorSlabs, Parks, Courtyards, Trees, Parcels, Metrics, MorphologyMetrics, and CaseManifest.
 
-`MorphologyMetrics` reports plan area density (`lambda_p`), open-space and park ratios, directional gross frontal area density (`lambda_f`), and building-height statistics including mean, standard deviation, median, and 90th percentile. `CaseManifest` is schema-versioned JSON containing a deterministic SHA-256 case ID, the effective generation parameters, random seed, plug-in version, street-network selection, model units, geometry counts, development metrics, and morphology metrics. Its schema is published at [`docs/case-manifest.schema.json`](docs/case-manifest.schema.json).
+`MorphologyMetrics` reports the urban form indices of Xu et al. (2017) and related work — built coverage ratio (BCR, identical to plan area density `lambda_p`), plan-area-weighted building height (BHt), built density ratio (BDR, built volume / site area), road area ratio (RaR), vertical area ratio (VR, facade area / site area), mean ground-level sky view factor (SVF, ray-cast over a grid of open-ground points against the masses and terrain), and aspect ratio (AR, mean height / `StreetWidth`) — alongside open-space and park ratios, directional gross frontal area density (`lambda_f`), and building-height statistics including mean, standard deviation, median, and 90th percentile. `CaseManifest` is schema-versioned JSON containing a deterministic SHA-256 case ID, the effective generation parameters, random seed, plug-in version, street-network selection, model units, geometry counts, development metrics, and morphology metrics. Its schema is published at [`docs/case-manifest.schema.json`](docs/case-manifest.schema.json).
 
 Each Building Type Config component exposes: floor range, corner radius, minimum footprint area, setback range, and building depth range. Feed any combination of configs into the generator — each block picks one at random.
 

docs/case-manifest.schema.json

@@ -7,7 +7,7 @@
   "required": ["schema", "schemaVersion", "caseId", "generator", "parameters", "geometry", "development", "morphology"],
   "properties": {
     "schema": { "type": "string", "format": "uri" },
-    "schemaVersion": { "const": "1.1.0" },
+    "schemaVersion": { "const": "1.2.0" },
     "caseId": { "type": "string", "pattern": "^[0-9a-f]{64}$" },
     "generator": {
       "type": "object",

src/Mycelium/Components/MassingGeneratorComponent.cs

@@ -82,7 +82,7 @@ protected override void RegisterOutputParams(GH_OutputParamManager pManager)
             pManager.AddCurveParameter("Parcels", "Parc", "Building parcel boundaries", GH_ParamAccess.list);
             pManager.AddTextParameter("Metrics", "Met", "Area and unit metrics", GH_ParamAccess.item);
             pManager.AddTextParameter("MorphologyMetrics", "Morph",
-                "Urban morphology indicators including lambda_p, directional lambda_f, open/park ratios, and height statistics",
+                "Urban morphology indicators including BCR (lambda_p), BDR, RaR, VR, SVF, AR, directional lambda_f, open/park ratios, and height statistics",
                 GH_ParamAccess.item);
             pManager.AddTextParameter("CaseManifest", "JSON",
                 "Versioned JSON case manifest containing inputs, provenance, counts, development metrics, and morphology metrics",
@@ -233,7 +233,8 @@ protected override void SolveInstance(IGH_DataAccess DA)
             // 5. Development and morphology metrics
             var development = CalculateDevelopmentMetrics(boundary, totalGFA);
             string metrics = BuildMetrics(development, allParcels.Count, masses.Count, parks.Count, trees.Count);
-            var morphology = MorphologyMetrics.Calculate(boundary, footprints, masses, parks, analysisDirection);
+            var morphology = MorphologyMetrics.Calculate(boundary, footprints, masses, parks, allParcels,
+                streetWidth, terrain, analysisDirection);
 
             // 6. Reproducible parameter/provenance sidecar. Geometry can be exported independently.
             var manifest = new CaseManifest

src/Mycelium/Core/CaseManifest.cs

@@ -17,7 +17,7 @@ namespace Mycelium.Core
     /// </summary>
     public sealed class CaseManifest
     {
-        public const string CurrentSchemaVersion = "1.1.0";
+        public const string CurrentSchemaVersion = "1.2.0";
 
         public string Schema { get; set; } = "https://github.com/MyceliumGH-Dev/Mycelium/blob/dev/docs/case-manifest.schema.json";
         public string SchemaVersion { get; set; } = CurrentSchemaVersion;

src/Mycelium/Core/MorphologyMetrics.cs

@@ -32,23 +32,40 @@ public sealed class MorphologyMetricsResult
         public double AnalysisDirectionY { get; set; }
         public int BuildingFootprintCount { get; set; }
         public int MassCount { get; set; }
+        public double BuiltVolume { get; set; }
+        public double BuiltDensityRatio { get; set; }
+        public double RoadArea { get; set; }
+        public double RoadAreaRatio { get; set; }
+        public double FacadeArea { get; set; }
+        public double VerticalAreaRatio { get; set; }
+        public double SkyViewFactor { get; set; }
+        public int SkyViewFactorSampleCount { get; set; }
+        public double AspectRatio { get; set; }
 
         public string ToDisplayString()
         {
             var text = new StringBuilder();
             text.AppendLine("--- Morphology Metrics ---");
-            text.AppendLine($"Plan Area Density (λp): {PlanAreaDensity:F3}");
+            text.AppendLine($"Built Coverage Ratio (BCR, λp): {PlanAreaDensity:F3}");
+            text.AppendLine($"Built Density Ratio (BDR): {BuiltDensityRatio:F3}");
+            text.AppendLine($"Road Area Ratio (RaR): {RoadAreaRatio:F3}");
+            text.AppendLine($"Vertical Area Ratio (VR): {VerticalAreaRatio:F3}");
+            text.AppendLine($"Sky View Factor (SVF, ground mean, {SkyViewFactorSampleCount} points): {SkyViewFactor:F3}");
+            text.AppendLine($"Aspect Ratio (AR, H/W): {AspectRatio:F3}");
             text.AppendLine($"Open Space Ratio: {OpenSpaceRatio:F3}");
             text.AppendLine($"Park Area Ratio: {ParkAreaRatio:F3}");
             text.AppendLine($"Frontal Area Density (λf, projected, unshielded): {FrontalAreaDensity:F3}");
             text.AppendLine($"Analysis Direction: ({AnalysisDirectionX:F3}, {AnalysisDirectionY:F3})");
             text.AppendLine($"Mean Height: {MeanHeight:F2}");
             text.AppendLine($"Height Std. Dev. (σH): {HeightStandardDeviation:F2}");
-            text.AppendLine($"Plan-Area-Weighted Mean Height: {PlanAreaWeightedMeanHeight:F2}");
+            text.AppendLine($"Building Height (BHt, plan-area-weighted mean): {PlanAreaWeightedMeanHeight:F2}");
             text.AppendLine($"Plan-Area-Weighted σH: {PlanAreaWeightedHeightStandardDeviation:F2}");
             text.AppendLine($"Height Min / Median / P90 / Max: {MinimumHeight:F2} / {MedianHeight:F2} / {HeightP90:F2} / {MaximumHeight:F2}");
             text.AppendLine($"Building Plan Area: {BuildingPlanArea:F2}");
             text.AppendLine($"Park Area: {ParkArea:F2}");
+            text.AppendLine($"Road Area: {RoadArea:F2}");
+            text.AppendLine($"Facade Area: {FacadeArea:F2}");
+            text.AppendLine($"Built Volume: {BuiltVolume:F2}");
             text.Append($"Gross Frontal Area: {GrossFrontalArea:F2}");
             return text.ToString();
         }
@@ -66,8 +83,21 @@ internal sealed class HeightStatistics
             public double P90 { get; set; }
         }
 
+        /// <summary>Target number of ground sample points for the sky view factor.</summary>
+        internal const int SkyViewSampleTarget = 400;
+
+        /// <summary>Rays cast per ground sample point for the sky view factor.</summary>
+        internal const int SkyViewRayCount = 144;
+
+        /// <param name="allParcels">
+        /// Every planar parcel the subdivision produced, park parcels included. Streets are the
+        /// site left over once these are removed, so road area is site area minus their sum.
+        /// </param>
+        /// <param name="streetWidth">Canyon width used for the aspect ratio.</param>
+        /// <param name="terrain">Optional terrain; sky view samples sit on it and it occludes.</param>
         public static MorphologyMetricsResult Calculate(Curve boundary, IReadOnlyList<Curve> footprints,
-            IReadOnlyList<Brep> masses, IReadOnlyList<Curve> parks, Vector3d analysisDirection)
+            IReadOnlyList<Brep> masses, IReadOnlyList<Curve> parks, IReadOnlyList<Curve> allParcels,
+            double streetWidth, Brep terrain, Vector3d analysisDirection)
         {
             double siteArea = GeometryHelpers.GetCurveArea(boundary);
             double planArea = GeometryHelpers.GetRegionArea(footprints);
@@ -84,6 +114,12 @@ public static MorphologyMetricsResult Calculate(Curve boundary, IReadOnlyList<Cu
             var weightedStatistics = CalculateWeightedHeightStatistics(weightedHeights);
             double grossFrontalArea = CalculateGrossFrontalArea(masses, direction);
 
+            double builtVolume = weightedHeights.Sum(entry => entry.PlanArea * entry.Height);
+            double roadArea = Math.Max(0.0, siteArea - SumCurveAreas(allParcels));
+            double facadeArea = CalculateFacadeArea(masses, weightedHeights);
+            double skyViewFactor = CalculateSkyViewFactor(boundary, siteArea, masses, terrain,
+                out int skyViewSamples);
+
             return new MorphologyMetricsResult
             {
                 SiteArea = siteArea,
@@ -105,7 +141,17 @@ public static MorphologyMetricsResult Calculate(Curve boundary, IReadOnlyList<Cu
                 AnalysisDirectionX = direction.X,
                 AnalysisDirectionY = direction.Y,
                 BuildingFootprintCount = footprints?.Count ?? 0,
-                MassCount = masses?.Count ?? 0
+                MassCount = masses?.Count ?? 0,
+                BuiltVolume = builtVolume,
+                BuiltDensityRatio = SafeRatio(builtVolume, siteArea),
+                RoadArea = roadArea,
+                RoadAreaRatio = SafeRatio(roadArea, siteArea),
+                FacadeArea = facadeArea,
+                VerticalAreaRatio = SafeRatio(facadeArea, siteArea),
+                SkyViewFactor = skyViewFactor,
+                SkyViewFactorSampleCount = skyViewSamples,
+                // Mean of H / W over buildings; W is the one street width the generator uses.
+                AspectRatio = SafeRatio(heightStatistics.Mean, streetWidth)
             };
         }
 
@@ -252,6 +298,161 @@ private static double CalculateGrossFrontalArea(IReadOnlyList<Brep> masses, Vect
             return total;
         }
 
+        /// <summary>
+        /// Facade (wall) area summed over masses: perimeter times height. Recovered as total
+        /// surface area minus roof and base, which is exact for the vertical prisms this
+        /// generator extrudes and counts courtyard walls as facade.
+        /// </summary>
+        /// <remarks>
+        /// On terrain the base of a mass is extended below ground as a foundation, so the
+        /// buried part of the wall is included.
+        /// </remarks>
+        private static double CalculateFacadeArea(IReadOnlyList<Brep> masses,
+            IReadOnlyList<WeightedHeight> weightedHeights)
+        {
+            if (masses == null)
+                return 0.0;
+
+            double surfaceArea = 0.0;
+            foreach (var mass in masses)
+            {
+                if (mass == null)
+                    continue;
+
+                using (var area = AreaMassProperties.Compute(mass))
+                {
+                    if (area != null && !double.IsNaN(area.Area) && !double.IsInfinity(area.Area))
+                        surfaceArea += area.Area;
+                }
+            }
+
+            return FacadeAreaFromSurface(surfaceArea, weightedHeights.Sum(entry => entry.PlanArea));
+        }
+
+        internal static double FacadeAreaFromSurface(double surfaceArea, double planArea)
+        {
+            return Math.Max(0.0, surfaceArea - 2.0 * Math.Max(0.0, planArea));
+        }
+
+        /// <summary>
+        /// Mean sky view factor over open ground inside the boundary. Sample points lie on a
+        /// regular grid; each casts a fixed cosine-weighted ray set against the masses (and the
+        /// terrain, if any), and its SVF is the unobstructed fraction. Points under a mass are
+        /// skipped, so the mean covers streets, parks, courtyards and setbacks.
+        /// </summary>
+        /// <remarks>
+        /// This is a geometric SVF for a horizontal plane at ground level. Trees are not
+        /// occluders, and buildings outside the boundary do not exist, so points near the edge
+        /// read higher than they would in a continuous city.
+        /// </remarks>
+        private static double CalculateSkyViewFactor(Curve boundary, double siteArea,
+            IReadOnlyList<Brep> masses, Brep terrain, out int sampleCount)
+        {
+            sampleCount = 0;
+            if (boundary == null || siteArea <= 0.0)
+                return 0.0;
+
+            var massMesh = MeshBreps(masses);
+            var terrainMesh = terrain != null ? MeshBreps(new[] { terrain }) : null;
+            var rays = CosineWeightedHemisphere(SkyViewRayCount);
+
+            var box = boundary.GetBoundingBox(true);
+            var plane = new Plane(new Point3d(0.0, 0.0, box.Min.Z), Vector3d.ZAxis);
+            double spacing = Math.Sqrt(siteArea / SkyViewSampleTarget);
+            double lift = spacing * 1e-3;
+            double skyTop = box.Max.Z;
+            if (massMesh != null)
+                skyTop = Math.Max(skyTop, massMesh.GetBoundingBox(false).Max.Z);
+            if (terrainMesh != null)
+                skyTop = Math.Max(skyTop, terrainMesh.GetBoundingBox(false).Max.Z);
+
+            double total = 0.0;
+            for (double x = box.Min.X + spacing / 2.0; x < box.Max.X; x += spacing)
+            {
+                for (double y = box.Min.Y + spacing / 2.0; y < box.Max.Y; y += spacing)
+                {
+                    var point = new Point3d(x, y, box.Min.Z);
+                    if (boundary.Contains(point, plane, 1e-6) != PointContainment.Inside)
+                        continue;
+
+                    if (terrainMesh != null)
+                    {
+                        var down = new Ray3d(new Point3d(x, y, skyTop + spacing), -Vector3d.ZAxis);
+                        double hit = Rhino.Geometry.Intersect.Intersection.MeshRay(terrainMesh, down);
+                        if (hit < 0.0)
+                            continue;
+                        point = down.PointAt(hit);
+                    }
+
+                    point.Z += lift;
+
+                    // A point under a mass sees a roof straight up; it is not open ground.
+                    if (massMesh != null &&
+                        Rhino.Geometry.Intersect.Intersection.MeshRay(massMesh, new Ray3d(point, Vector3d.ZAxis)) >= 0.0)
+                        continue;
+
+                    int open = 0;
+                    foreach (var direction in rays)
+                    {
+                        var ray = new Ray3d(point, direction);
+                        bool blocked =
+                            (massMesh != null && Rhino.Geometry.Intersect.Intersection.MeshRay(massMesh, ray) >= 0.0) ||
+                            (terrainMesh != null && Rhino.Geometry.Intersect.Intersection.MeshRay(terrainMesh, ray) >= 0.0);
+                        if (!blocked)
+                            open++;
+                    }
+
+                    total += (double)open / rays.Length;
+                    sampleCount++;
+                }
+            }
+
+            return sampleCount > 0 ? total / sampleCount : 0.0;
+        }
+
+        /// <summary>
+        /// Deterministic cosine-weighted directions over the upper hemisphere (Fibonacci spiral
+        /// on the unit disk, lifted to the sphere). With this density every ray carries equal
+        /// weight, so the sky view factor is simply the unobstructed fraction.
+        /// </summary>
+        internal static Vector3d[] CosineWeightedHemisphere(int count)
+        {
+            var directions = new Vector3d[Math.Max(0, count)];
+            double goldenAngle = Math.PI * (3.0 - Math.Sqrt(5.0));
+            for (int i = 0; i < directions.Length; i++)
+            {
+                double u = (i + 0.5) / directions.Length;
+                double radius = Math.Sqrt(u);
+                double angle = i * goldenAngle;
+                directions[i] = new Vector3d(radius * Math.Cos(angle), radius * Math.Sin(angle),
+                    Math.Sqrt(1.0 - u));
+            }
+
+            return directions;
+        }
+
+        private static Mesh MeshBreps(IEnumerable<Brep> breps)
+        {
+            if (breps == null)
+                return null;
+
+            var joined = new Mesh();
+            foreach (var brep in breps)
+            {
+                if (brep == null)
+                    continue;
+
+                var pieces = Mesh.CreateFromBrep(brep, MeshingParameters.FastRenderMesh);
+                if (pieces == null)
+                    continue;
+
+                foreach (var piece in pieces)
+                    joined.Append(piece);
+            }
+
+            return joined.Faces.Count > 0 ? joined : null;
+        }
+
         private static double SumCurveAreas(IReadOnlyList<Curve> curves)
         {
             if (curves == null)

tests/Mycelium.Core.Tests/MorphologyMetricsTests.cs

@@ -51,6 +51,36 @@ public void WeightedHeightStatistics_FallBackWhenNoPlanAreaIsAvailable()
             Assert.Equal(15.0, MorphologyMetrics.CalculateWeightedHeightStatistics(entries).Mean, 8);
         }
 
+        [Fact]
+        public void FacadeArea_IsSurfaceMinusRoofAndBase()
+        {
+            // 10 x 10 x 20 prism: surface 1000, plan 100, walls 4 * 10 * 20 = 800.
+            Assert.Equal(800.0, MorphologyMetrics.FacadeAreaFromSurface(1000.0, 100.0), 8);
+            Assert.Equal(0.0, MorphologyMetrics.FacadeAreaFromSurface(100.0, 100.0), 8);
+        }
+
+        [Fact]
+        public void CosineWeightedHemisphere_IsUnitUpwardAndCosineDistributed()
+        {
+            var rays = MorphologyMetrics.CosineWeightedHemisphere(MorphologyMetrics.SkyViewRayCount);
+
+            Assert.Equal(MorphologyMetrics.SkyViewRayCount, rays.Length);
+            double meanZ = 0.0, meanX = 0.0, meanY = 0.0;
+            foreach (var ray in rays)
+            {
+                Assert.True(ray.Z > 0.0);
+                Assert.Equal(1.0, ray.X * ray.X + ray.Y * ray.Y + ray.Z * ray.Z, 10);
+                meanX += ray.X / rays.Length;
+                meanY += ray.Y / rays.Length;
+                meanZ += ray.Z / rays.Length;
+            }
+
+            // E[cos(theta)] under a cosine-weighted density is 2/3; the set has no azimuthal bias.
+            Assert.Equal(2.0 / 3.0, meanZ, 2);
+            Assert.Equal(0.0, meanX, 2);
+            Assert.Equal(0.0, meanY, 2);
+        }
+
         [Fact]
         public void CaseManifest_UsesStableCamelCaseSchema()
         {