feat(tdigest): implement tdigest.cdf command

src/commands/cmd_tdigest.cc

@@ -556,6 +556,47 @@ class CommandTDigestTrimmedMean : public Commander {
   double high_cut_quantile_;
 };
 
+class CommandTDigestCDF : public Commander {
+  Status Parse(const std::vector<std::string> &args) override {
+    if (args.size() == 2) return {Status::RedisParseErr, errWrongNumOfArguments};
+    key_name_ = args[1];
+    inputs_.reserve(args.size() - 2);
+    for (size_t i = 2; i < args.size(); i++) {
+      auto value = ParseFloat(args[i]);
+      if (!value) {
+        return {Status::RedisParseErr, errValueIsNotFloat};
+      }
+      if (std::isnan(*value)) {
+        return {Status::RedisParseErr, errValueIsNotFloat};
+      }
+      inputs_.push_back(*value);
+    }
+    return Status::OK();
+  }
+
+  Status Execute(engine::Context &ctx, Server *srv, Connection *conn, std::string *output) override {
+    TDigest tdigest(srv->storage, conn->GetNamespace());
+    TDigestCDFResult result;
+    auto s = tdigest.CDF(ctx, key_name_, inputs_, &result);
+    if (!s.ok()) {
+      if (s.IsNotFound()) {
+        return {Status::RedisExecErr, errKeyNotFound};
+      }
+      return {Status::RedisExecErr, s.ToString()};
+    }
+
+    output->append(redis::MultiLen(result.cdf_values.size()));
+    for (auto const value : result.cdf_values) {
+      output->append(conn->Double(value));
+    }
+    return Status::OK();
+  }
+
+ private:
+  std::string key_name_;
+  std::vector<double> inputs_;
+};
+
 std::vector<CommandKeyRange> GetMergeKeyRange(const std::vector<std::string> &args) {
   auto numkeys = ParseInt<int>(args[2], 10).ValueOr(0);
   return {{1, 1, 1}, {3, 2 + numkeys, 1}};
@@ -573,5 +614,6 @@ REDIS_REGISTER_COMMANDS(TDigest, MakeCmdAttr<CommandTDigestCreate>("tdigest.crea
                         MakeCmdAttr<CommandTDigestQuantile>("tdigest.quantile", -3, "read-only", 1, 1, 1),
                         MakeCmdAttr<CommandTDigestTrimmedMean>("tdigest.trimmed_mean", 4, "read-only", 1, 1, 1),
                         MakeCmdAttr<CommandTDigestReset>("tdigest.reset", 2, "write", 1, 1, 1),
-                        MakeCmdAttr<CommandTDigestMerge>("tdigest.merge", -4, "write", GetMergeKeyRange));
+                        MakeCmdAttr<CommandTDigestMerge>("tdigest.merge", -4, "write", GetMergeKeyRange),
+                        MakeCmdAttr<CommandTDigestCDF>("tdigest.cdf", -3, "read-only", 1, 1, 1));
 }  // namespace redis

src/types/redis_tdigest.cc

@@ -28,12 +28,14 @@
 #include <rocksdb/status.h>
 
 #include <algorithm>
+#include <cstdint>
 #include <iterator>
 #include <limits>
 #include <memory>
 #include <range/v3/algorithm/minmax.hpp>
 #include <range/v3/range/conversion.hpp>
 #include <range/v3/view/join.hpp>
+#include <range/v3/view/map.hpp>
 #include <range/v3/view/transform.hpp>
 #include <vector>
 
@@ -570,6 +572,41 @@ rocksdb::Status TDigest::Merge(engine::Context& ctx, const Slice& dest_digest,
   return storage_->Write(ctx, storage_->DefaultWriteOptions(), batch->GetWriteBatch());
 }
 
+rocksdb::Status TDigest::CDF(engine::Context& ctx, const Slice& digest_name, const std::vector<double>& inputs,
+                             TDigestCDFResult* result) {
+  auto ns_key = AppendNamespacePrefix(digest_name);
+  TDigestMetadata metadata;
+  {
+    LockGuard guard(storage_->GetLockManager(), ns_key);
+
+    if (auto status = getMetaDataByNsKey(ctx, ns_key, &metadata); !status.ok()) {
+      return status;
+    }
+
+    if (metadata.total_observations == 0) {
+      result->cdf_values = std::vector<double>(inputs.size(), std::numeric_limits<double>::quiet_NaN());
+      return rocksdb::Status::OK();
+    }
+
+    if (auto status = mergeNodes(ctx, ns_key, &metadata); !status.ok()) {
+      return status;
+    }
+  }
+
+  std::vector<Centroid> centroids;
+  if (auto status = dumpCentroids(ctx, ns_key, metadata, &centroids); !status.ok()) {
+    return status;
+  }
+
+  auto dump_centroids = DummyCentroids<false>(metadata, centroids);
+  if (auto status = TDigestCDF(centroids, dump_centroids.Min(), dump_centroids.Max(), dump_centroids.TotalWeight(),
+                               inputs, &result->cdf_values);
+      !status.IsOK()) {
+    return rocksdb::Status::InvalidArgument(status.Msg());
+  }
+  return rocksdb::Status::OK();
+}
+
 rocksdb::Status TDigest::GetMetaData(engine::Context& context, const Slice& digest_name, TDigestMetadata* metadata) {
   auto ns_key = AppendNamespacePrefix(digest_name);
   return Database::GetMetadata(context, {kRedisTDigest}, ns_key, metadata);

src/types/redis_tdigest.h

@@ -49,6 +49,10 @@ struct TDigestMergeOptions {
   bool override_flag = false;
 };
 
+struct TDigestCDFResult {
+  std::vector<double> cdf_values;
+};
+
 struct TDigestQuantitleResult {
   std::optional<std::vector<double>> quantiles;
 };
@@ -93,6 +97,9 @@ class TDigest : public SubKeyScanner {
                               double high_cut_quantile, TDigestTrimmedMeanResult* result);
   rocksdb::Status GetMetaData(engine::Context& context, const Slice& digest_name, TDigestMetadata* metadata);
 
+  rocksdb::Status CDF(engine::Context& ctx, const Slice& digest_name, const std::vector<double>& inputs,
+                      TDigestCDFResult* result);
+
  private:
   enum class SegmentType : uint8_t { kBuffer = 0, kCentroids = 1, kGuardFlag = 0xFF };
 

src/types/tdigest.h

@@ -22,10 +22,11 @@
 
 #include <fmt/format.h>
 
+#include <algorithm>
+#include <cmath>
+#include <iterator>
 #include <limits>
 #include <map>
-#include <numeric>
-#include <variant>
 #include <vector>
 
 #include "common/status.h"
@@ -172,6 +173,131 @@ struct DoubleComparator {
   bool operator()(const double& a, const double& b) const { return DoubleCompare(a, b) == -1; }
 };
 
+inline double InterpolateRank(double value, double lower_value, double upper_value, double lower_rank,
+                              double upper_rank) {
+  if (DoubleEqual(lower_value, upper_value)) {
+    return (lower_rank + upper_rank) / 2;
+  }
+  return Lerp(lower_rank, upper_rank, (value - lower_value) / (upper_value - lower_value));
+}
+
+// Match RedisBloom t-digest-c CDF behavior: if min/max is outside the first/last centroid mean, the exact
+// boundary sample is treated as a singleton with weight 1. Its center rank is 0.5 at min and
+// total_weight - 0.5 at max; interpolation toward an inner centroid starts after the singleton, at rank 1 or
+// total_weight - 1.
+inline Status TDigestCDF(const std::vector<Centroid>& centroids, double min, double max, double total_weight,
+                         const std::vector<double>& inputs, std::vector<double>* result) {
+  if (centroids.empty() || total_weight <= 0) {
+    return Status{Status::InvalidArgument, "invalid or empty tdigest"};
+  }
+
+  std::map<double, std::vector<size_t>> sorted_unique_input_idx_map;
+  for (size_t i = 0; i < inputs.size(); ++i) {
+    sorted_unique_input_idx_map[inputs[i]].push_back(i);
+  }
+
+  std::vector<double> sorted_unique_inputs;
+  sorted_unique_inputs.reserve(sorted_unique_input_idx_map.size());
+  std::transform(sorted_unique_input_idx_map.cbegin(), sorted_unique_input_idx_map.cend(),
+                 std::back_inserter(sorted_unique_inputs), [](const auto& pair) { return pair.first; });
+
+  constexpr double kSingletonBoundaryWeight = 1.0;
+  constexpr double kHalfSingletonBoundaryWeight = kSingletonBoundaryWeight / 2;
+
+  struct CentroidGroup {
+    double mean;
+    double weight;
+    double center_rank;
+  };
+
+  std::vector<CentroidGroup> groups;
+  groups.reserve(centroids.size());
+  double cumulative_weight = 0;
+  for (size_t i = 0; i < centroids.size();) {
+    double group_weight = 0;
+    const auto mean = centroids[i].mean;
+    do {
+      group_weight += centroids[i].weight;
+      ++i;
+    } while (i < centroids.size() && DoubleEqual(mean, centroids[i].mean));
+
+    groups.push_back({
+        .mean = mean,
+        .weight = group_weight,
+        .center_rank = cumulative_weight + group_weight / 2,
+    });
+    cumulative_weight += group_weight;
+  }
+
+  std::vector<double> sorted_results;
+  sorted_results.reserve(sorted_unique_inputs.size());
+  size_t group_idx = 0;
+  for (const auto value : sorted_unique_inputs) {
+    if (value < min) {
+      sorted_results.push_back(0);
+      continue;
+    }
+    if (value > max) {
+      sorted_results.push_back(1);
+      continue;
+    }
+
+    if (value == min) {
+      auto rank = DoubleEqual(groups.front().mean, min) ? groups.front().center_rank : kHalfSingletonBoundaryWeight;
+      sorted_results.push_back(rank / total_weight);
+      continue;
+    }
+    if (value == max) {
+      auto rank = DoubleEqual(groups.back().mean, max) ? groups.back().center_rank
+                                                       : (total_weight - kHalfSingletonBoundaryWeight);
+      sorted_results.push_back(rank / total_weight);
+      continue;
+    }
+
+    while (group_idx < groups.size() && DoubleCompare(groups[group_idx].mean, value) < 0) {
+      ++group_idx;
+    }
+
+    double rank = 0;
+    if (group_idx == groups.size()) {
+      const auto& last = groups.back();
+      rank = InterpolateRank(value, last.mean, max, last.center_rank, total_weight - kSingletonBoundaryWeight);
+    } else if (DoubleEqual(value, groups[group_idx].mean)) {
+      rank = groups[group_idx].center_rank;
+    } else if (group_idx == 0) {
+      rank = InterpolateRank(value, min, groups.front().mean, kSingletonBoundaryWeight, groups.front().center_rank);
+    } else {
+      const auto& left = groups[group_idx - 1];
+      const auto& right = groups[group_idx];
+      auto weight_before_left = left.center_rank - left.weight / 2;
+      // When both adjacent centroids are singletons, include the left exact sample and exclude the right one.
+      if (left.weight == kSingletonBoundaryWeight && right.weight == kSingletonBoundaryWeight) {
+        rank = weight_before_left + kSingletonBoundaryWeight;
+      } else {
+        // Exclude singleton half-weights from interpolation because singleton centroids are exact samples.
+        double left_excluded_weight = left.weight == kSingletonBoundaryWeight ? kHalfSingletonBoundaryWeight : 0;
+        double right_excluded_weight = right.weight == kSingletonBoundaryWeight ? kHalfSingletonBoundaryWeight : 0;
+        double rank_span = (left.weight + right.weight) / 2 - left_excluded_weight - right_excluded_weight;
+        double base_rank = weight_before_left + left.weight / 2 + left_excluded_weight;
+        rank = base_rank + rank_span * (value - left.mean) / (right.mean - left.mean);
+      }
+    }
+
+    rank = std::clamp(rank, 0.0, total_weight);
+    sorted_results.push_back(rank / total_weight);
+  }
+
+  result->clear();
+  result->resize(inputs.size(), std::numeric_limits<double>::quiet_NaN());
+  for (size_t i = 0; i < sorted_unique_inputs.size(); ++i) {
+    for (auto idx : sorted_unique_input_idx_map[sorted_unique_inputs[i]]) {
+      (*result)[idx] = sorted_results[i];
+    }
+  }
+
+  return Status::OK();
+}
+
 template <bool Reverse, typename TD>
 inline Status TDigestByRank(TD&& td, const std::vector<int>& inputs, std::vector<double>* result) {
   result->clear();