Geometry Nodes: Attribute Statistic Node

This nodes evaluates a field on a geometry and outputs various
statistics about the entire data set, like min, max, or even
the standard deviation. It works for float and vector types currently,
though more types could be supported in the future.

 - All statistics are calculated element-wise for vectors.
 - "Product" was not added since the result could very easily overflow.
 - The "Size" output was not added since it isn't specific to an
   attribute and would fit better in another node.

The implementation shares work as much as possible when multiple
statistics are needed.

This node has been in development since the beginning of this year,
with additions from Johnny Matthews and Hans Goudey.

Differential Revision: https://developer.blender.org/D10202

release/scripts/startup/nodeitems_builtins.py

@@ -503,6 +503,7 @@ geometry_node_categories = [
         NodeItem("GeometryNodeAttributeRemove", poll=geometry_nodes_fields_legacy_poll),
 
         NodeItem("GeometryNodeAttributeCapture", poll=geometry_nodes_fields_poll),
+        NodeItem("GeometryNodeAttributeStatistic", poll=geometry_nodes_fields_poll),
     ]),
     GeometryNodeCategory("GEO_COLOR", "Color", items=[
         NodeItem("ShaderNodeMixRGB"),

source/blender/blenkernel/BKE_node.h

@@ -1494,6 +1494,7 @@ int ntreeTexExecTree(struct bNodeTree *ntree,
 #define GEO_NODE_MATERIAL_SELECTION 1081
 #define GEO_NODE_MATERIAL_ASSIGN 1082
 #define GEO_NODE_REALIZE_INSTANCES 1083
+#define GEO_NODE_ATTRIBUTE_STATISTIC 1084
 
 /** \} */
 

source/blender/blenkernel/intern/node.cc

@@ -5175,6 +5175,7 @@ static void registerGeometryNodes()
   register_node_type_geo_attribute_randomize();
   register_node_type_geo_attribute_remove();
   register_node_type_geo_attribute_separate_xyz();
+  register_node_type_geo_attribute_statistic();
   register_node_type_geo_attribute_transfer();
   register_node_type_geo_attribute_vector_math();
   register_node_type_geo_attribute_vector_rotate();

source/blender/blenlib/BLI_array.hh

@@ -276,6 +276,21 @@ class Array {
     initialized_fill_n(data_, size_, value);
   }
 
+  /**
+   * Return a reference to the first element in the array.
+   * This invokes undefined behavior when the array is empty.
+   */
+  const T &first() const
+  {
+    BLI_assert(size_ > 0);
+    return *data_;
+  }
+  T &first()
+  {
+    BLI_assert(size_ > 0);
+    return *data_;
+  }
+
   /**
    * Return a reference to the last element in the array.
    * This invokes undefined behavior when the array is empty.

source/blender/makesrna/intern/rna_nodetree.c

@@ -2168,6 +2168,17 @@ static const EnumPropertyItem *rna_GeometryNodeAttributeFill_type_itemf(bContext
   return itemf_function_check(rna_enum_attribute_type_items, attribute_fill_type_supported);
 }
 
+static bool attribute_statistic_type_supported(const EnumPropertyItem *item)
+{
+  return ELEM(item->value, CD_PROP_FLOAT, CD_PROP_FLOAT3);
+}
+static const EnumPropertyItem *rna_GeometryNodeAttributeStatistic_type_itemf(
+    bContext *UNUSED(C), PointerRNA *UNUSED(ptr), PropertyRNA *UNUSED(prop), bool *r_free)
+{
+  *r_free = true;
+  return itemf_function_check(rna_enum_attribute_type_items, attribute_statistic_type_supported);
+}
+
 /**
  * This bit of ugly code makes sure the float / attribute option shows up instead of
  * vector / attribute if the node uses an operation that uses a float for input B or C.
@@ -9219,6 +9230,29 @@ static void def_geo_attribute_convert(StructRNA *srna)
   RNA_def_property_update(prop, NC_NODE | NA_EDITED, "rna_Node_update");
 }
 
+static void def_geo_attribute_statistic(StructRNA *srna)
+{
+  PropertyRNA *prop;
+
+  prop = RNA_def_property(srna, "data_type", PROP_ENUM, PROP_NONE);
+  RNA_def_property_enum_sdna(prop, NULL, "custom1");
+  RNA_def_property_enum_items(prop, rna_enum_attribute_type_items);
+  RNA_def_property_enum_funcs(prop, NULL, NULL, "rna_GeometryNodeAttributeStatistic_type_itemf");
+  RNA_def_property_enum_default(prop, CD_PROP_FLOAT);
+  RNA_def_property_ui_text(
+      prop,
+      "Data Type",
+      "The data type the attribute is converted to before calculating the results");
+  RNA_def_property_update(prop, NC_NODE | NA_EDITED, "rna_GeometryNode_socket_update");
+
+  prop = RNA_def_property(srna, "domain", PROP_ENUM, PROP_NONE);
+  RNA_def_property_enum_sdna(prop, NULL, "custom2");
+  RNA_def_property_enum_items(prop, rna_enum_attribute_domain_items);
+  RNA_def_property_enum_default(prop, ATTR_DOMAIN_POINT);
+  RNA_def_property_ui_text(prop, "Domain", "Which domain to read the data from");
+  RNA_def_property_update(prop, NC_NODE | NA_EDITED, "rna_Node_update");
+}
+
 static void def_geo_attribute_math(StructRNA *srna)
 {
   PropertyRNA *prop;

source/blender/nodes/CMakeLists.txt

@@ -161,6 +161,7 @@ set(SRC
   geometry/nodes/node_geo_attribute_remove.cc
   geometry/nodes/node_geo_attribute_sample_texture.cc
   geometry/nodes/node_geo_attribute_separate_xyz.cc
+  geometry/nodes/node_geo_attribute_statistic.cc
   geometry/nodes/node_geo_attribute_transfer.cc
   geometry/nodes/node_geo_attribute_vector_math.cc
   geometry/nodes/node_geo_attribute_vector_rotate.cc

source/blender/nodes/NOD_geometry.h

@@ -48,6 +48,7 @@ void register_node_type_geo_attribute_proximity(void);
 void register_node_type_geo_attribute_randomize(void);
 void register_node_type_geo_attribute_remove(void);
 void register_node_type_geo_attribute_separate_xyz(void);
+void register_node_type_geo_attribute_statistic(void);
 void register_node_type_geo_attribute_transfer(void);
 void register_node_type_geo_attribute_vector_math(void);
 void register_node_type_geo_attribute_vector_rotate(void);

source/blender/nodes/NOD_static_types.h

@@ -307,6 +307,7 @@ DefNode(GeometryNode, GEO_NODE_LEGACY_SELECT_BY_MATERIAL, 0, "LEGACY_SELECT_BY_M
 DefNode(GeometryNode, GEO_NODE_ATTRIBUTE_CAPTURE, def_geo_attribute_capture, "ATTRIBUTE_CAPTURE", AttributeCapture, "Attribute Capture", "")
 DefNode(GeometryNode, GEO_NODE_ATTRIBUTE_REMOVE, 0, "ATTRIBUTE_REMOVE", AttributeRemove, "Attribute Remove", "")
 DefNode(GeometryNode, GEO_NODE_ATTRIBUTE_VECTOR_ROTATE, def_geo_attribute_vector_rotate, "LEGACY_ATTRIBUTE_VECTOR_ROTATE", LegacyAttributeVectorRotate, "Attribute Vector Rotate", "")
+DefNode(GeometryNode, GEO_NODE_ATTRIBUTE_STATISTIC, def_geo_attribute_statistic, "ATTRIBUTE_STATISTIC", AttributeStatistic, "Attribute Statistic", "")
 DefNode(GeometryNode, GEO_NODE_BOOLEAN, def_geo_boolean, "BOOLEAN", Boolean, "Boolean", "")
 DefNode(GeometryNode, GEO_NODE_BOUNDING_BOX, 0, "BOUNDING_BOX", BoundBox, "Bounding Box", "")
 DefNode(GeometryNode, GEO_NODE_COLLECTION_INFO, def_geo_collection_info, "COLLECTION_INFO", CollectionInfo, "Collection Info", "")

source/blender/nodes/geometry/nodes/node_geo_attribute_statistic.cc

@@ -0,0 +1,378 @@
+/*
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ */
+
+#include <algorithm>
+#include <numeric>
+
+#include "UI_interface.h"
+#include "UI_resources.h"
+
+#include "BLI_math_base_safe.h"
+
+#include "node_geometry_util.hh"
+
+namespace blender::nodes {
+
+static void geo_node_attribute_statistic_declare(NodeDeclarationBuilder &b)
+{
+  b.add_input<decl::Geometry>("Geometry");
+  b.add_input<decl::Float>("Attribute").hide_value();
+  b.add_input<decl::Vector>("Attribute", "Attribute_001").hide_value();
+
+  b.add_output<decl::Float>("Mean");
+  b.add_output<decl::Float>("Median");
+  b.add_output<decl::Float>("Sum");
+  b.add_output<decl::Float>("Min");
+  b.add_output<decl::Float>("Max");
+  b.add_output<decl::Float>("Range");
+  b.add_output<decl::Float>("Standard Deviation");
+  b.add_output<decl::Float>("Variance");
+
+  b.add_output<decl::Vector>("Mean", "Mean_001");
+  b.add_output<decl::Vector>("Median", "Median_001");
+  b.add_output<decl::Vector>("Sum", "Sum_001");
+  b.add_output<decl::Vector>("Min", "Min_001");
+  b.add_output<decl::Vector>("Max", "Max_001");
+  b.add_output<decl::Vector>("Range", "Range_001");
+  b.add_output<decl::Vector>("Standard Deviation", "Standard Deviation_001");
+  b.add_output<decl::Vector>("Variance", "Variance_001");
+}
+
+static void geo_node_attribute_statistic_layout(uiLayout *layout,
+                                                bContext *UNUSED(C),
+                                                PointerRNA *ptr)
+{
+  uiItemR(layout, ptr, "data_type", 0, "", ICON_NONE);
+  uiItemR(layout, ptr, "domain", 0, "", ICON_NONE);
+}
+
+static void geo_node_attribute_statistic_init(bNodeTree *UNUSED(tree), bNode *node)
+{
+  node->custom1 = CD_PROP_FLOAT;
+  node->custom2 = ATTR_DOMAIN_POINT;
+}
+
+static void geo_node_attribute_statistic_update(bNodeTree *UNUSED(ntree), bNode *node)
+{
+  bNodeSocket *socket_geo = (bNodeSocket *)node->inputs.first;
+  bNodeSocket *socket_float_attr = socket_geo->next;
+  bNodeSocket *socket_float3_attr = socket_float_attr->next;
+
+  bNodeSocket *socket_float_mean = (bNodeSocket *)node->outputs.first;
+  bNodeSocket *socket_float_median = socket_float_mean->next;
+  bNodeSocket *socket_float_sum = socket_float_median->next;
+  bNodeSocket *socket_float_min = socket_float_sum->next;
+  bNodeSocket *socket_float_max = socket_float_min->next;
+  bNodeSocket *socket_float_range = socket_float_max->next;
+  bNodeSocket *socket_float_std = socket_float_range->next;
+  bNodeSocket *socket_float_variance = socket_float_std->next;
+
+  bNodeSocket *socket_vector_mean = socket_float_variance->next;
+  bNodeSocket *socket_vector_median = socket_vector_mean->next;
+  bNodeSocket *socket_vector_sum = socket_vector_median->next;
+  bNodeSocket *socket_vector_min = socket_vector_sum->next;
+  bNodeSocket *socket_vector_max = socket_vector_min->next;
+  bNodeSocket *socket_vector_range = socket_vector_max->next;
+  bNodeSocket *socket_vector_std = socket_vector_range->next;
+  bNodeSocket *socket_vector_variance = socket_vector_std->next;
+
+  const CustomDataType data_type = static_cast<CustomDataType>(node->custom1);
+
+  nodeSetSocketAvailability(socket_float_attr, data_type == CD_PROP_FLOAT);
+  nodeSetSocketAvailability(socket_float_mean, data_type == CD_PROP_FLOAT);
+  nodeSetSocketAvailability(socket_float_median, data_type == CD_PROP_FLOAT);
+  nodeSetSocketAvailability(socket_float_sum, data_type == CD_PROP_FLOAT);
+  nodeSetSocketAvailability(socket_float_min, data_type == CD_PROP_FLOAT);
+  nodeSetSocketAvailability(socket_float_max, data_type == CD_PROP_FLOAT);
+  nodeSetSocketAvailability(socket_float_range, data_type == CD_PROP_FLOAT);
+  nodeSetSocketAvailability(socket_float_std, data_type == CD_PROP_FLOAT);
+  nodeSetSocketAvailability(socket_float_variance, data_type == CD_PROP_FLOAT);
+
+  nodeSetSocketAvailability(socket_float3_attr, data_type == CD_PROP_FLOAT3);
+  nodeSetSocketAvailability(socket_vector_mean, data_type == CD_PROP_FLOAT3);
+  nodeSetSocketAvailability(socket_vector_median, data_type == CD_PROP_FLOAT3);
+  nodeSetSocketAvailability(socket_vector_sum, data_type == CD_PROP_FLOAT3);
+  nodeSetSocketAvailability(socket_vector_min, data_type == CD_PROP_FLOAT3);
+  nodeSetSocketAvailability(socket_vector_max, data_type == CD_PROP_FLOAT3);
+  nodeSetSocketAvailability(socket_vector_range, data_type == CD_PROP_FLOAT3);
+  nodeSetSocketAvailability(socket_vector_std, data_type == CD_PROP_FLOAT3);
+  nodeSetSocketAvailability(socket_vector_variance, data_type == CD_PROP_FLOAT3);
+}
+
+template<typename T> static T compute_sum(const Span<T> data)
+{
+  return std::accumulate(data.begin(), data.end(), T());
+}
+
+static float compute_variance(const Span<float> data, const float mean)
+{
+  if (data.size() <= 1) {
+    return 0.0f;
+  }
+
+  float sum_of_squared_differences = std::accumulate(
+      data.begin(), data.end(), 0.0f, [mean](float accumulator, float value) {
+        float difference = mean - value;
+        return accumulator + difference * difference;
+      });
+
+  return sum_of_squared_differences / (data.size() - 1);
+}
+
+static float median_of_sorted_span(const Span<float> data)
+{
+  if (data.is_empty()) {
+    return 0.0f;
+  }
+
+  const float median = data[data.size() / 2];
+
+  /* For spans of even length, the median is the average of the middle two elements. */
+  if (data.size() % 2 == 0) {
+    return (median + data[data.size() / 2 - 1]) * 0.5f;
+  }
+  return median;
+}
+static void set_empty(CustomDataType data_type, GeoNodeExecParams &params)
+{
+  if (data_type == CD_PROP_FLOAT) {
+    params.set_output("Mean", 0.0f);
+    params.set_output("Median", 0.0f);
+    params.set_output("Sum", 0.0f);
+    params.set_output("Min", 0.0f);
+    params.set_output("Max", 0.0f);
+    params.set_output("Range", 0.0f);
+    params.set_output("Standard Deviation", 0.0f);
+    params.set_output("Variance", 0.0f);
+  }
+  else if (data_type == CD_PROP_FLOAT3) {
+    params.set_output("Mean_001", float3{0.0f, 0.0f, 0.0f});
+    params.set_output("Median_001", float3{0.0f, 0.0f, 0.0f});
+    params.set_output("Sum_001", float3{0.0f, 0.0f, 0.0f});
+    params.set_output("Min_001", float3{0.0f, 0.0f, 0.0f});
+    params.set_output("Max_001", float3{0.0f, 0.0f, 0.0f});
+    params.set_output("Range_001", float3{0.0f, 0.0f, 0.0f});
+    params.set_output("Standard Deviation_001", float3{0.0f, 0.0f, 0.0f});
+    params.set_output("Variance_001", float3{0.0f, 0.0f, 0.0f});
+  }
+}
+
+static void geo_node_attribute_statistic_exec(GeoNodeExecParams params)
+{
+  GeometrySet geometry_set = params.get_input<GeometrySet>("Geometry");
+
+  const bNode &node = params.node();
+  const CustomDataType data_type = static_cast<CustomDataType>(node.custom1);
+  const AttributeDomain domain = static_cast<AttributeDomain>(node.custom2);
+
+  int64_t total_size = 0;
+  Vector<const GeometryComponent *> components = geometry_set.get_components_for_read();
+
+  for (const GeometryComponent *component : components) {
+    if (component->attribute_domain_supported(domain)) {
+      total_size += component->attribute_domain_size(domain);
+    }
+  }
+  if (total_size == 0) {
+    set_empty(data_type, params);
+    return;
+  }
+
+  switch (data_type) {
+    case CD_PROP_FLOAT: {
+      const Field<float> input_field = params.get_input<Field<float>>("Attribute");
+      Array<float> data = Array<float>(total_size);
+      int offset = 0;
+      for (const GeometryComponent *component : components) {
+        if (component->attribute_domain_supported(domain)) {
+          GeometryComponentFieldContext field_context{*component, domain};
+          const int domain_size = component->attribute_domain_size(domain);
+          fn::FieldEvaluator data_evaluator{field_context, domain_size};
+          MutableSpan<float> component_result = data.as_mutable_span().slice(offset, domain_size);
+          data_evaluator.add_with_destination(input_field, component_result);
+          data_evaluator.evaluate();
+          offset += domain_size;
+        }
+      }
+
+      float mean = 0.0f;
+      float median = 0.0f;
+      float sum = 0.0f;
+      float min = 0.0f;
+      float max = 0.0f;
+      float range = 0.0f;
+      float standard_deviation = 0.0f;
+      float variance = 0.0f;
+      const bool sort_required = params.output_is_required("Min") ||
+                                 params.output_is_required("Max") ||
+                                 params.output_is_required("Range") ||
+                                 params.output_is_required("Median");
+      const bool sum_required = params.output_is_required("Sum") ||
+                                params.output_is_required("Mean");
+      const bool variance_required = params.output_is_required("Standard Deviation") ||
+                                     params.output_is_required("Variance");
+
+      if (total_size != 0) {
+        if (sort_required) {
+          std::sort(data.begin(), data.end());
+          median = median_of_sorted_span(data);
+
+          min = data.first();
+          max = data.last();
+          range = max - min;
+        }
+        if (sum_required || variance_required) {
+          sum = compute_sum<float>(data);
+          mean = sum / total_size;
+
+          if (variance_required) {
+            variance = compute_variance(data, mean);
+            standard_deviation = std::sqrt(variance);
+          }
+        }
+      }
+
+      if (sum_required) {
+        params.set_output("Sum", sum);
+        params.set_output("Mean", mean);
+      }
+      if (sort_required) {
+        params.set_output("Min", min);
+        params.set_output("Max", max);
+        params.set_output("Range", range);
+        params.set_output("Median", median);
+      }
+      if (variance_required) {
+        params.set_output("Standard Deviation", standard_deviation);
+        params.set_output("Variance", variance);
+      }
+      break;
+    }
+    case CD_PROP_FLOAT3: {
+      const Field<float3> input_field = params.get_input<Field<float3>>("Attribute_001");
+
+      Array<float3> data = Array<float3>(total_size);
+      int offset = 0;
+      for (const GeometryComponent *component : components) {
+        if (component->attribute_domain_supported(domain)) {
+          GeometryComponentFieldContext field_context{*component, domain};
+          const int domain_size = component->attribute_domain_size(domain);
+          fn::FieldEvaluator data_evaluator{field_context, domain_size};
+          MutableSpan<float3> component_result = data.as_mutable_span().slice(offset, domain_size);
+          data_evaluator.add_with_destination(input_field, component_result);
+          data_evaluator.evaluate();
+          offset += domain_size;
+        }
+      }
+
+      float3 median{0};
+      float3 min{0};
+      float3 max{0};
+      float3 range{0};
+      float3 sum{0};
+      float3 mean{0};
+      float3 variance{0};
+      float3 standard_deviation{0};
+      const bool sort_required = params.output_is_required("Min_001") ||
+                                 params.output_is_required("Max_001") ||
+                                 params.output_is_required("Range_001") ||
+                                 params.output_is_required("Median_001");
+      const bool sum_required = params.output_is_required("Sum_001") ||
+                                params.output_is_required("Mean_001");
+      const bool variance_required = params.output_is_required("Standard Deviation_001") ||
+                                     params.output_is_required("Variance_001");
+
+      Array<float> data_x;
+      Array<float> data_y;
+      Array<float> data_z;
+      if (sort_required || variance_required) {
+        data_x.reinitialize(total_size);
+        data_y.reinitialize(total_size);
+        data_z.reinitialize(total_size);
+        for (const int i : data.index_range()) {
+          data_x[i] = data[i].x;
+          data_y[i] = data[i].y;
+          data_z[i] = data[i].z;
+        }
+      }
+
+      if (total_size != 0) {
+        if (sort_required) {
+          std::sort(data_x.begin(), data_x.end());
+          std::sort(data_y.begin(), data_y.end());
+          std::sort(data_z.begin(), data_z.end());
+
+          const float x_median = median_of_sorted_span(data_x);
+          const float y_median = median_of_sorted_span(data_y);
+          const float z_median = median_of_sorted_span(data_z);
+          median = float3(x_median, y_median, z_median);
+
+          min = float3(data_x.first(), data_y.first(), data_z.first());
+          max = float3(data_x.last(), data_y.last(), data_z.last());
+          range = max - min;
+        }
+        if (sum_required || variance_required) {
+          sum = compute_sum(data.as_span());
+          mean = sum / total_size;
+
+          if (variance_required) {
+            const float x_variance = compute_variance(data_x, mean.x);
+            const float y_variance = compute_variance(data_y, mean.y);
+            const float z_variance = compute_variance(data_z, mean.z);
+            variance = float3(x_variance, y_variance, z_variance);
+            standard_deviation = float3(
+                std::sqrt(variance.x), std::sqrt(variance.y), std::sqrt(variance.z));
+          }
+        }
+      }
+
+      if (sum_required) {
+        params.set_output("Sum_001", sum);
+        params.set_output("Mean_001", mean);
+      }
+      if (sort_required) {
+        params.set_output("Min_001", min);
+        params.set_output("Max_001", max);
+        params.set_output("Range_001", range);
+        params.set_output("Median_001", median);
+      }
+      if (variance_required) {
+        params.set_output("Standard Deviation_001", standard_deviation);
+        params.set_output("Variance_001", variance);
+      }
+      break;
+    }
+    default:
+      break;
+  }
+}
+
+}  // namespace blender::nodes
+
+void register_node_type_geo_attribute_statistic()
+{
+  static bNodeType ntype;
+
+  geo_node_type_base(
+      &ntype, GEO_NODE_ATTRIBUTE_STATISTIC, "Attribute Statistic", NODE_CLASS_ATTRIBUTE, 0);
+
+  ntype.declare = blender::nodes::geo_node_attribute_statistic_declare;
+  node_type_init(&ntype, blender::nodes::geo_node_attribute_statistic_init);
+  node_type_update(&ntype, blender::nodes::geo_node_attribute_statistic_update);
+  ntype.geometry_node_execute = blender::nodes::geo_node_attribute_statistic_exec;
+  ntype.draw_buttons = blender::nodes::geo_node_attribute_statistic_layout;
+  nodeRegisterType(&ntype);
+}