#include <bits/stdc++.h>
#include "artclass.h"
using namespace std;
constexpr double DThreshold = 200;
constexpr size_t FRatio = 2;
vector<vector<int64_t>> sobel = vector<vector<int64_t>>({vector<int64_t>({-1, 0, 1}),
vector<int64_t>({-2, 0, 2}),
vector<int64_t>({-1, 0, 1})});
void bit_reverse(vector<complex<double>> &y, size_t lgn)
{
vector<bool> swapped(y.size(), 0);
for (size_t i = 0; i < y.size(); i++)
{
if (swapped[i])
continue;
size_t j = 0;
for (size_t k = 0; k < lgn; k++)
j |= ((i & (1 << k)) >> k) << (lgn - k - 1);
swap(y[i], y[j]);
swapped[i] = swapped[j] = 1;
}
}
void fft(vector<complex<double>> &y)
{
size_t const n = y.size(), lgn = __lg(n);
bit_reverse(y, lgn);
for (size_t s = 1; s <= lgn; s++)
{
size_t const m = 1 << s;
complex<double> const omega_m = polar(1.0, 2 * M_PI / m);
for (size_t i = 0; i < n; i += m)
{
complex<double> omega = 1.0;
for (size_t k = i; k < i + m / 2; k++)
{
complex<double> const u = y[k], v = omega * y[k + m / 2];
y[k] = u + v;
y[k + m / 2] = u - v;
omega *= omega_m;
}
}
}
}
template <typename T>
void rotate_c(vector<vector<T>> &y)
{
vector<vector<T>> z(y[0].size(), vector<T>(y.size()));
for (size_t i = 0; i < y.size(); i++)
for (size_t j = 0; j < y[0].size(); j++)
z[j][y.size() - i - 1] = y[i][j];
swap(y, z);
}
template <typename T>
size_t convolve(vector<vector<T>> const &y, vector<vector<T>> const &kernel)
{
size_t counter = 0;
for (size_t i = kernel.size() / 2; i < y.size() - kernel.size() / 2; i++)
for (size_t j = kernel[0].size() / 2; j < y[i].size() - kernel[0].size() / 2; j++)
{
T val = 0;
for (size_t k = 0; k < kernel.size(); k++)
for (size_t h = 0; h < kernel[k].size(); h++)
val += y[i - kernel.size() / 2 + k][j - kernel.size() / 2 + h] * kernel[k][h];
if (abs(val) > DThreshold)
counter++;
}
return counter;
}
size_t get_high_deriv(int H, int W, int arr[500][500])
{
vector<vector<int64_t>> image(H, vector<int64_t>(W));
for (size_t i = 0; i < H; i++)
for (size_t j = 0; j < W; j++)
image[i][j] = arr[i][j];
size_t ans = convolve(image, sobel);
rotate_c(sobel);
ans += convolve(image, sobel);
rotate_c(sobel);
rotate_c(sobel);
rotate_c(sobel);
return ans;
}
double get_high_freq(int H, int W, int arr[500][500])
{
size_t n1 = 1 << (32 - __builtin_clz(H)), n2 = 1 << (32 - __builtin_clz(W));
vector<vector<complex<double>>> y(n1, vector<complex<double>>(n2));
for (size_t i = 0; i < H; i++)
for (size_t j = 0; j < W; j++)
y[i][j] = arr[i][j];
for (size_t i = 0; i < n1; i++)
fft(y[i]);
rotate_c(y);
for (size_t i = 0; i < n2; i++)
fft(y[i]);
double ans = 0.0;
for (size_t i = n2 - n2 / FRatio; i < n2; i++)
for (size_t j = n1 - n1 / FRatio; j < n1; j++)
ans += abs(y[i][j]);
return ans;
}
int style(int H, int W, int R[500][500], int G[500][500], int B[500][500])
{
double h = (double)(get_high_deriv(H, W, R) + get_high_deriv(H, W, G) +
get_high_deriv(H, W, B)) /
(double)(H * W);
double f = (get_high_freq(H, W, R) + get_high_freq(H, W, G) + get_high_freq(H, W, B)) / (H * W);
// cout << "h : " << h << ", "
// << "f : " << f << '\n';
if (h < 0.15)
{
if (h < 0.06)
return f < 2200 ? 4 : 1;
return (f < 2500 || h < 0.01) ? 4 : 2;
}
return (f > 10000 || f > 7000 && h > 0.7) ? 3 : (f / h > 25000 ? 2 : 1);
}
