Mirror/aicodix___dsp
1
0
Fork 0
mirror of https://github.com/aicodix/dsp.git synced 2026-04-27 14:30:36 +00:00
Code Issues Projects Releases Packages Wiki Activity

use quick select to get the median

Browse source
This commit is contained in:
Ahmet Inan 2024-03-09 15:07:48 +01:00
commit 46403a7b0a
2 changed files with 11 additions and 20 deletions
Download patch file Download diff file Expand all files Collapse all files

23
repeated_median.hh
View file

@ -6,7 +6,7 @@ Copyright 2021 Ahmet Inan <inan@aicodix.de>
#pragma once #pragma once
#include <algorithm> #include "quick_select.hh"
namespace DSP { namespace DSP {
@ -26,21 +26,17 @@ public:
for (int j = 0; j < LEN; ++j) for (int j = 0; j < LEN; ++j)
if (x[j] != x[i]) if (x[j] != x[i])
inner_[count++] = (y[j] - y[i]) / (x[j] - x[i]); inner_[count++] = (y[j] - y[i]) / (x[j] - x[i]);
std::nth_element(inner_, inner_+count/2, inner_+count); outer_[i] = quick_select(inner_, count/2, count);
outer_[i] = inner_[count/2];
} }
std::nth_element(outer_, outer_+LEN/2, outer_+LEN); slope_ = quick_select(outer_, LEN/2, LEN);
slope_ = outer_[LEN/2];
for (int i = 0; i < LEN; ++i) { for (int i = 0; i < LEN; ++i) {
int count = 0; int count = 0;
for (int j = 0; j < LEN; ++j) for (int j = 0; j < LEN; ++j)
if (x[j] != x[i]) if (x[j] != x[i])
inner_[count++] = (x[j]*y[i] - x[i]*y[j]) / (x[j] - x[i]); inner_[count++] = (x[j]*y[i] - x[i]*y[j]) / (x[j] - x[i]);
std::nth_element(inner_, inner_+count/2, inner_+count); outer_[i] = quick_select(inner_, count/2, count);
outer_[i] = inner_[count/2];
} }
std::nth_element(outer_, outer_+LEN/2, outer_+LEN); yint_ = quick_select(outer_, LEN/2, LEN);
yint_ = outer_[LEN/2];
xint_ = - yint_ / slope_; xint_ = - yint_ / slope_;
} }
TYPE xint() TYPE xint()
@ -77,15 +73,12 @@ public:
for (int j = 0; j < LEN; ++j) for (int j = 0; j < LEN; ++j)
if (x[j] != x[i]) if (x[j] != x[i])
inner_[count++] = (y[j] - y[i]) / (x[j] - x[i]); inner_[count++] = (y[j] - y[i]) / (x[j] - x[i]);
std::nth_element(inner_, inner_+count/2, inner_+count); outer_[i] = quick_select(inner_, count/2, count);
outer_[i] = inner_[count/2];
} }
std::nth_element(outer_, outer_+LEN/2, outer_+LEN); slope_ = quick_select(outer_, LEN/2, LEN);
slope_ = outer_[LEN/2];
for (int i = 0; i < LEN; ++i) for (int i = 0; i < LEN; ++i)
outer_[i] = y[i] - slope_ * x[i]; outer_[i] = y[i] - slope_ * x[i];
std::nth_element(outer_, outer_+LEN/2, outer_+LEN); yint_ = quick_select(outer_, LEN/2, LEN);
yint_ = outer_[LEN/2];
xint_ = - yint_ / slope_; xint_ = - yint_ / slope_;
} }
TYPE xint() TYPE xint()

8
theil_sen.hh
View file

@ -6,7 +6,7 @@ Copyright 2021 Ahmet Inan <inan@aicodix.de>
#pragma once #pragma once
#include <algorithm> #include "quick_select.hh"
namespace DSP { namespace DSP {
@ -25,13 +25,11 @@ public:
for (int j = i+1; count < size_ && j < LEN; ++j) for (int j = i+1; count < size_ && j < LEN; ++j)
if (x[j] != x[i]) if (x[j] != x[i])
temp_[count++] = (y[j] - y[i]) / (x[j] - x[i]); temp_[count++] = (y[j] - y[i]) / (x[j] - x[i]);
std::nth_element(temp_, temp_+count/2, temp_+count); slope_ = quick_select(temp_, count/2, count);
slope_ = temp_[count/2];
count = 0; count = 0;
for (int i = 0; count < size_ && i < LEN; ++i) for (int i = 0; count < size_ && i < LEN; ++i)
temp_[count++] = y[i] - slope_ * x[i]; temp_[count++] = y[i] - slope_ * x[i];
std::nth_element(temp_, temp_+count/2, temp_+count); yint_ = quick_select(temp_, count/2, count);
yint_ = temp_[count/2];
xint_ = - yint_ / slope_; xint_ = - yint_ / slope_;
} }
TYPE xint() TYPE xint()