//+------------------------------------------------------------------+
//|                              BenchSubrangeAB_AllInOne.mq5         |
//| Foco: SOLO el caso interesante -> m_reconstruir=true con subrango |
//| real (fs != total_cols), comparando:                              |
//|   CPU  (patron validado: out=X + overwrite del rango)             |
//|   GPU-A (mandar TODO + offsets, leer completo, sin reconstruccion |
//|          manual porque el resto nunca se toco)                    |
//|   GPU-B (mandar SOLO el rango, reconstruir con out=X + overwrite) |
//|                                                                    |
//| El caso fs==total_cols y el caso !reconstruir NO se benchmarkean  |
//| aca a proposito: ya sabemos que van por la ruta trivial (0 copias |
//| intermedias, sin ambiguedad de estrategia).                       |
//+------------------------------------------------------------------+
#property script_show_inputs

#ifndef SCALERBYLEO_CL_DEVICE
#define SCALERBYLEO_CL_DEVICE (CL_USE_GPU_ONLY)
#endif

//+------------------------------------------------------------------+
static const string MinMaxKernelSrcTemplate =
   "#define MATRIX_IDX(row, col, cols) ((row) * (cols) + (col)) \r\n"
   "\r\n"
   "__kernel void max_min_scaler_transform_matrix( \r\n"
   "   __global const real_t* min_values, \r\n"
   "   __global const real_t* diff_inv, \r\n"
   "   __global real_t* matrix, \r\n"
   "   const uint cols, \r\n"
   "   const uint istart \r\n"
   ") \r\n"
   "{ \r\n"
   "   const uint row = get_global_id(0); \r\n"
   "   const uint col = get_global_id(1); \r\n"
   "   const uint matrix_idx = MATRIX_IDX(row, istart + col, cols); \r\n"
   "   matrix[matrix_idx] = (matrix[matrix_idx] - min_values[col]) * diff_inv[col]; \r\n"
   "} \r\n"
   "\r\n"
   "__kernel void max_min_scaler_transform_matrix_exact( \r\n"
   "   __global const real_t* min_values, \r\n"
   "   __global const real_t* diff_inv, \r\n"
   "   __global real_t* matrix, \r\n"
   "   const uint cols \r\n"
   ") \r\n"
   "{ \r\n"
   "   const uint row = get_global_id(0); \r\n"
   "   const uint col = get_global_id(1); \r\n"
   "   const uint matrix_idx = MATRIX_IDX(row, col, cols); \r\n"
   "   matrix[matrix_idx] = (matrix[matrix_idx] - min_values[col]) * diff_inv[col]; \r\n"
   "} \r\n";

//+------------------------------------------------------------------+
class CCLCtxShared
 {
private:
   static int s_ref_count;
public:
   static int s_context;

   static bool ContextInit(void)
    {
      if(s_context != INVALID_HANDLE) { s_ref_count++; return true; }
      ::ResetLastError();
      if((s_context = CLContextCreate(SCALERBYLEO_CL_DEVICE)) == INVALID_HANDLE)
         return false;
      s_ref_count = 1;
      return true;
    }

   static void ContextRelease(void)
    {
      if(s_context == INVALID_HANDLE) return;
      if(--s_ref_count <= 0)
       {
         CLContextFree(s_context);
         s_context = INVALID_HANDLE;
         s_ref_count = 0;
       }
    }
 };
int CCLCtxShared::s_context   = INVALID_HANDLE;
int CCLCtxShared::s_ref_count = 0;

//+------------------------------------------------------------------+
template<typename TFlt>
class CCLProgramManager
 {
public:
   static int  s_program;
   static int  s_kernel_full;
   static int  s_kernel_exact;
   static bool s_ready;

   static bool Init(void)
    {
      if(s_ready) return true;
      if(!CCLCtxShared::ContextInit()) return false;

      const string real_t_def = (sizeof(TFlt) == sizeof(float))
                                 ? "#define real_t float \r\n"
                                 : "#define real_t double \r\n";
      const string src = real_t_def + MinMaxKernelSrcTemplate;

      string build_log;
      ::ResetLastError();
      if((s_program = CLProgramCreate(CCLCtxShared::s_context, src, build_log)) == INVALID_HANDLE)
       {
         Print("CLProgramCreate failed: ", build_log, " err=", GetLastError());
         CCLCtxShared::ContextRelease();
         return false;
       }
      if((s_kernel_full = CLKernelCreate(s_program, "max_min_scaler_transform_matrix")) == INVALID_HANDLE)
       {
         Print("CLKernelCreate(full) failed, err=", GetLastError());
         return false;
       }
      if((s_kernel_exact = CLKernelCreate(s_program, "max_min_scaler_transform_matrix_exact")) == INVALID_HANDLE)
       {
         Print("CLKernelCreate(exact) failed, err=", GetLastError());
         return false;
       }
      s_ready = true;
      return true;
    }

   static void Deinit(void)
    {
      if(!s_ready) return;
      if(s_kernel_full != INVALID_HANDLE)  { CLKernelFree(s_kernel_full);  s_kernel_full  = INVALID_HANDLE; }
      if(s_kernel_exact != INVALID_HANDLE) { CLKernelFree(s_kernel_exact); s_kernel_exact = INVALID_HANDLE; }
      if(s_program != INVALID_HANDLE)      { CLProgramFree(s_program);    s_program      = INVALID_HANDLE; }
      CCLCtxShared::ContextRelease();
      s_ready = false;
    }
 };
template<typename TFlt> int  CCLProgramManager::s_program     = INVALID_HANDLE;
template<typename TFlt> int  CCLProgramManager::s_kernel_full  = INVALID_HANDLE;
template<typename TFlt> int  CCLProgramManager::s_kernel_exact = INVALID_HANDLE;
template<typename TFlt> bool CCLProgramManager::s_ready        = false;

//+------------------------------------------------------------------+
//| MinMax scaler: ruta trivial + dos estrategias para subrango real |
//+------------------------------------------------------------------+
template<typename TFlt>
class MaxMinScalerCL
 {
private:
   bool         m_loaded_scaler;
   ulong        m_i_start;
   vector<TFlt> m_min_vals;
   vector<TFlt> m_diff_inv;
   int          m_buffer_h_min_vals;
   int          m_buffer_h_diff_inv;

   bool EnsureAuxBuffers(void)
    {
      const uint bytes = uint(m_min_vals.Size() * sizeof(TFlt));
      if(m_buffer_h_min_vals == INVALID_HANDLE)
         if((m_buffer_h_min_vals = CLBufferCreate(CCLCtxShared::s_context, bytes, CL_MEM_READ_ONLY)) == INVALID_HANDLE)
            return false;
      if(!CLBufferWrite(m_buffer_h_min_vals, 0, m_min_vals)) return false;

      if(m_buffer_h_diff_inv == INVALID_HANDLE)
         if((m_buffer_h_diff_inv = CLBufferCreate(CCLCtxShared::s_context, bytes, CL_MEM_READ_ONLY)) == INVALID_HANDLE)
            return false;
      if(!CLBufferWrite(m_buffer_h_diff_inv, 0, m_diff_inv)) return false;

      return true;
    }

public:
   MaxMinScalerCL(void) : m_loaded_scaler(false), m_i_start(0),
                           m_buffer_h_min_vals(INVALID_HANDLE), m_buffer_h_diff_inv(INVALID_HANDLE) {}

   ~MaxMinScalerCL(void)
    {
      if(m_buffer_h_min_vals != INVALID_HANDLE) CLBufferFree(m_buffer_h_min_vals);
      if(m_buffer_h_diff_inv != INVALID_HANDLE) CLBufferFree(m_buffer_h_diff_inv);
    }

   bool Fit(const matrix<TFlt> &X, const ulong i_start, const ulong fs)
    {
      if(!CCLProgramManager<TFlt>::Init()) return false;
      m_i_start = i_start;

      m_min_vals.Resize(fs);
      m_diff_inv.Resize(fs);

      for(ulong c = 0; c < fs; c++)
       {
         TFlt mn = X[0][i_start + c];
         TFlt mx = X[0][i_start + c];
         for(ulong r = 1; r < X.Rows(); r++)
          {
            const TFlt v = X[r][i_start + c];
            if(v < mn) mn = v;
            if(v > mx) mx = v;
          }
         const TFlt diff = (mx - mn);
         m_min_vals[c] = mn;
         m_diff_inv[c] = (diff == 0) ? TFlt(0) : TFlt(1.0 / double(diff));
       }

      m_loaded_scaler = true;
      return EnsureAuxBuffers();
    }

   //--- SUBRANGO ESTRATEGIA A: mandar TODO + offsets, leer completo directo (0 copias CPU extra)
   bool SubrangeStrategyA(matrix<TFlt> &X, matrix<TFlt> &out)
    {
      if(!m_loaded_scaler) return false;

      const ulong rows       = X.Rows();
      const ulong total_cols = X.Cols();
      const ulong fs         = m_min_vals.Size();

      int matrix_h = INVALID_HANDLE;
      const uint bytes = uint(sizeof(TFlt) * total_cols * rows);
      matrix_h = CLBufferCreate(CCLCtxShared::s_context, bytes, CL_MEM_READ_WRITE);
      if(matrix_h == INVALID_HANDLE) return false;

      bool ok = CLBufferWrite(matrix_h, 0, X);

      if(ok)
         ok = CLSetKernelArgMem(CCLProgramManager<TFlt>::s_kernel_full, 0, m_buffer_h_min_vals) &&
              CLSetKernelArgMem(CCLProgramManager<TFlt>::s_kernel_full, 1, m_buffer_h_diff_inv) &&
              CLSetKernelArgMem(CCLProgramManager<TFlt>::s_kernel_full, 2, matrix_h)            &&
              CLSetKernelArg(CCLProgramManager<TFlt>::s_kernel_full, 3, uint(fs))                &&
              CLSetKernelArg(CCLProgramManager<TFlt>::s_kernel_full, 4, uint(m_i_start));

      if(ok)
       {
         const uint work_off[2] = {0, 0};
         uint work_size[2] = { (uint)rows, (uint)fs };
         ok = CLExecute(CCLProgramManager<TFlt>::s_kernel_full, 2, work_off, work_size);
       }

      if(ok)
       {
         out.Resize(rows, total_cols);
         ok = CLBufferRead(matrix_h, 0, out);
       }

      CLBufferFree(matrix_h);
      return ok;
    }

   //--- SUBRANGO ESTRATEGIA B: mandar SOLO el rango, reconstruir con out=X + overwrite
   bool SubrangeStrategyB(matrix<TFlt> &X, matrix<TFlt> &out)
    {
      if(!m_loaded_scaler) return false;

      const ulong rows = X.Rows();
      const ulong fs   = m_min_vals.Size();

      matrix<TFlt> narrow;
      narrow.Resize(rows, fs);
      for(ulong c = 0; c < fs; c++)
         for(ulong r = 0; r < rows; r++)
            narrow[r][c] = X[r][m_i_start + c];

      int matrix_h = INVALID_HANDLE;
      const uint bytes = uint(sizeof(TFlt) * fs * rows);
      matrix_h = CLBufferCreate(CCLCtxShared::s_context, bytes, CL_MEM_READ_WRITE);
      bool ok = (matrix_h != INVALID_HANDLE) && CLBufferWrite(matrix_h, 0, narrow);

      if(ok)
         ok = CLSetKernelArgMem(CCLProgramManager<TFlt>::s_kernel_exact, 0, m_buffer_h_min_vals) &&
              CLSetKernelArgMem(CCLProgramManager<TFlt>::s_kernel_exact, 1, m_buffer_h_diff_inv) &&
              CLSetKernelArgMem(CCLProgramManager<TFlt>::s_kernel_exact, 2, matrix_h)            &&
              CLSetKernelArg(CCLProgramManager<TFlt>::s_kernel_exact, 3, uint(fs));

      if(ok)
       {
         const uint work_off[2] = {0, 0};
         uint work_size[2] = { (uint)rows, (uint)fs };
         ok = CLExecute(CCLProgramManager<TFlt>::s_kernel_exact, 2, work_off, work_size);
       }

      if(!ok) { CLBufferFree(matrix_h); return false; }

      matrix<TFlt> res;
      res.Resize(rows, fs);
      ok = CLBufferRead(matrix_h, 0, res);
      CLBufferFree(matrix_h);
      if(!ok) return false;

      out = X;
      for(ulong c = 0; c < fs; c++)
         for(ulong r = 0; r < rows; r++)
            out[r][m_i_start + c] = res[r][c];

      return true;
    }
 };

//+------------------------------------------------------------------+
//| Version CPU pura (patron validado: out=X + overwrite del rango)  |
//+------------------------------------------------------------------+
template<typename TFlt>
class MaxMinScalerCPU
 {
private:
   ulong        m_i_start;
   vector<TFlt> m_min_vals;
   vector<TFlt> m_diff_inv;

public:
   MaxMinScalerCPU(void) : m_i_start(0) {}

   bool Fit(const matrix<TFlt> &X, const ulong i_start, const ulong fs)
    {
      m_i_start = i_start;
      m_min_vals.Resize(fs);
      m_diff_inv.Resize(fs);

      for(ulong c = 0; c < fs; c++)
       {
         TFlt mn = X[0][i_start + c];
         TFlt mx = X[0][i_start + c];
         for(ulong r = 1; r < X.Rows(); r++)
          {
            const TFlt v = X[r][i_start + c];
            if(v < mn) mn = v;
            if(v > mx) mx = v;
          }
         const TFlt diff = (mx - mn);
         m_min_vals[c] = mn;
         m_diff_inv[c] = (diff == 0) ? TFlt(0) : TFlt(1.0 / double(diff));
       }
      return true;
    }

   bool Transform(const matrix<TFlt> &X, matrix<TFlt> &out)
    {
      const ulong rows = X.Rows();
      const ulong fs    = m_min_vals.Size();

      out = X;
      for(ulong c = 0; c < fs; c++)
         for(ulong r = 0; r < rows; r++)
            out[r][m_i_start + c] = TFlt((double(X[r][m_i_start + c]) - double(m_min_vals[c])) * double(m_diff_inv[c]));

      return true;
    }
 };

//+------------------------------------------------------------------+
//| Config del barrido: SOLO subrango real, fracciones chicas a 90%  |
//+------------------------------------------------------------------+
int    g_sizes[]     = {500, 2000, 5000, 9600};
double g_fs_fracs[]   = {0.01, 0.05, 0.10, 0.25, 0.50, 0.75, 0.90};
int    g_repeats      = 5;

//+------------------------------------------------------------------+
ulong RunOnce_CPU(matrix<float> &X, const ulong i_start, const ulong fs)
 {
  MaxMinScalerCPU<float> sc;
  sc.Fit(X, i_start, fs);
  matrix<float> out;
  const ulong t0 = GetMicrosecondCount();
  sc.Transform(X, out);
  return GetMicrosecondCount() - t0;
 }

//+------------------------------------------------------------------+
ulong RunOnce_A(matrix<float> &X, const ulong i_start, const ulong fs)
 {
  MaxMinScalerCL<float> sc;
  sc.Fit(X, i_start, fs);
  matrix<float> out;
  const ulong t0 = GetMicrosecondCount();
  sc.SubrangeStrategyA(X, out);
  return GetMicrosecondCount() - t0;
 }

//+------------------------------------------------------------------+
ulong RunOnce_B(matrix<float> &X, const ulong i_start, const ulong fs)
 {
  MaxMinScalerCL<float> sc;
  sc.Fit(X, i_start, fs);
  matrix<float> out;
  const ulong t0 = GetMicrosecondCount();
  sc.SubrangeStrategyB(X, out);
  return GetMicrosecondCount() - t0;
 }

//+------------------------------------------------------------------+
ulong Median(ulong &arr[])
 {
  const int n = ArraySize(arr);
  ArraySort(arr);
  return arr[n / 2];
 }

//+------------------------------------------------------------------+
void OnStart()
 {
  Print("=== Bench SUBRANGO REAL (reconstruir=true, fs != total_cols): CPU vs GPU-A vs GPU-B ===");
  PrintFormat("%-6s %-6s %-6s %-12s %-12s %-12s %-10s %-10s",
              "size", "fs", "frac", "CPU(us)", "GPU_A(us)", "GPU_B(us)", "A_vs_CPU", "B_vs_CPU");

  for(int si = 0; si < ArraySize(g_sizes); si++)
   {
    const int size = g_sizes[si];

    matrix<float> X;
    X.Resize(size, size);
    X.Random(0.0, 20.0);

    for(int fi = 0; fi < ArraySize(g_fs_fracs); fi++)
     {
      ulong fs = (ulong)MathMax(1, MathRound(size * g_fs_fracs[fi]));
      if(fs >= (ulong)size) fs = size - 1; // forzamos subrango real, nunca fs==total_cols aca
      const ulong i_start = 0;

      ulong cpu_t[], a_t[], b_t[];
      ArrayResize(cpu_t, g_repeats);
      ArrayResize(a_t,   g_repeats);
      ArrayResize(b_t,   g_repeats);

      for(int rep = 0; rep < g_repeats; rep++)
       {
        cpu_t[rep] = RunOnce_CPU(X, i_start, fs);
        a_t[rep]   = RunOnce_A(X, i_start, fs);
        b_t[rep]   = RunOnce_B(X, i_start, fs);
       }

      const ulong cpu_med = Median(cpu_t);
      const ulong a_med   = Median(a_t);
      const ulong b_med   = Median(b_t);

      const double a_ratio = (a_med > 0) ? (double)cpu_med / (double)a_med : 0.0;
      const double b_ratio = (b_med > 0) ? (double)cpu_med / (double)b_med : 0.0;

      PrintFormat("%-6d %-6d %-6.2f %-12I64u %-12I64u %-12I64u %-10.2f %-10.2f",
                  size, fs, g_fs_fracs[fi], cpu_med, a_med, b_med, a_ratio, b_ratio);
     }
   }

  Print("=== Fin ===");
  Print("Ratio > 1.0 = GPU mas rapida que CPU. Comparar tambien GPU_A vs GPU_B directamente en la tabla.");

  CCLProgramManager<float>::Deinit();
 }
//+------------------------------------------------------------------+