new: add quick_select

src/quick_select.c

@@ -0,0 +1,55 @@
+#include <string.h>
+#include <stdlib.h>
+#include "quick_select.h"
+
+static inline void swap_double(double *a, double *b) {
+  double t = *a;
+  *a = *b;
+  *b = t;
+}
+
+// Lomuto partition
+static size_t partition(double *arr, size_t lo, size_t hi) {
+  double pivot = arr[hi];
+  size_t i = lo;
+  for (size_t j = lo; j < hi; ++j) {
+    if (arr[j] < pivot) {
+      swap_double(&arr[i], &arr[j]);
+      ++i;
+    }
+  }
+  swap_double(&arr[i], &arr[hi]);
+  return i;
+}
+
+// 选择第 k 小的元素 (0-based)，平均 O(n)，最坏 O(n^2)
+static double select_k(double *arr, size_t n, size_t k) {
+    size_t lo = 0;
+    size_t hi = n - 1;
+
+  for (;;) {
+    size_t p = partition(arr, lo, hi);
+    if (p == k) {
+      return arr[p];
+    } else if (k < p) {
+      hi = p - 1;
+    } else {
+      lo = p + 1;
+    }
+  }
+  return arr[k];
+}
+
+double percentile(const double *arr, size_t n, double perc) {
+  if (n == 0) {
+    return 0.0;
+  }
+
+  double *arr_copy = malloc(sizeof(double)*n);
+  memcpy(arr_copy, arr, n);
+  // nearest-rank: ceil(0.9*n) - 1
+  size_t k = (size_t)(perc*n) - 1;
+  double ans = select_k(arr_copy, n, k);
+  free(arr_copy);
+  return ans;
+}

src/quick_select.h

@@ -0,0 +1,2 @@
+#include <stddef.h>
+double percentile(const double *arr, size_t n, double perc);

src/write_png.c

@@ -4,6 +4,7 @@
 #include <stdio.h>
 #include "write_png.h"
 #include "common.h"
+#include "quick_select.h"
 
 static inline double linear_to_srgb(double x){
   if (x <= 0.0) return 0.0;
@@ -14,13 +15,14 @@ static inline double linear_to_srgb(double x){
 
 int buffer_normalize_srgb(double *buffer, int W, int H) {
   double max = 0;
-  for (int j = 0; j < H; j++) {
-    for (int i = 0; i < W; i++) {
-      for (int c = 0; c < 3; c++) {
-        if (max < buffer[color_index(i, j, c)]) max = buffer[color_index(i, j, c)];
-      }
-    }
-  }
+  //for (int j = 0; j < H; j++) {
+  //  for (int i = 0; i < W; i++) {
+  //    for (int c = 0; c < 3; c++) {
+  //      if (max < buffer[color_index(i, j, c)]) max = buffer[color_index(i, j, c)];
+  //    }
+  //  }
+  //}
+  max = percentile(buffer, W*H*3, cutperc);
   printf("max in the buffer: %g\n", max);
   if (max == 0) return 1;
   for (int j = 0; j < H; j++) {
@@ -48,7 +50,7 @@ int write_png(char *filename, double *buffer, int W, int H) {
         if (v < 0) v = 0;
         if (v > 255) v = 255;
         img[color_index(i, j, c)] = (uint8_t) v;
-        printf("i: %d, j: %d, c:%d, buffer=%g, v=%d\n", i,j,c,buffer[color_index(i, j, c)],v);
+        //printf("i: %d, j: %d, c:%d, buffer=%g, v=%d\n", i,j,c,buffer[color_index(i, j, c)],v);
       }
     }
   }