#include <iostream>
#include <vector>
#include <cstdint>
#include <algorithm>

#define STB_IMAGE_IMPLEMENTATION
#include "stb_image.h"

int main() {
    int width, height, channels;
    uint8_t* img = stbi_load("output2g.png", &width, &height, &channels, 1);
    if (!img) {
        std::cerr << "Failed to load output2g.png\n";
        return 1;
    }
    
    std::cout << "Image dimensions: " << height << " rows x " << width << " columns\n";
    
    std::vector<int> row_weights(height, 0);
    std::vector<int> col_weights(width, 0);
    
    for (int r = 0; r < height; r++) {
        for (int c = 0; c < width; c++) {
            if (img[r * width + c] < 128) {
                row_weights[r]++;
                col_weights[c]++;
            }
        }
    }
    
    int min_row = width, max_row = 0;
    long long sum_row = 0;
    for (int w : row_weights) {
        if (w < min_row) min_row = w;
        if (w > max_row) max_row = w;
        sum_row += w;
    }
    std::cout << "Row weights: min=" << min_row << ", max=" << max_row << ", avg=" << (sum_row / (double)height) << "\n";
    
    int min_col = height, max_col = 0;
    long long sum_col = 0;
    for (int w : col_weights) {
        if (w < min_col) min_col = w;
        if (w > max_col) max_col = w;
        sum_col += w;
    }
    std::cout << "Col weights: min=" << min_col << ", max=" << max_col << ", avg=" << (sum_col / (double)width) << "\n";
    
    int block_h = 256, block_w = 256;
    int num_blocks_y = height / block_h;
    int num_blocks_x = width / block_w;
    
    std::cout << "Analyzing " << num_blocks_y << "x" << num_blocks_x << " blocks of size 256x256\n";
    
    bool is_perfectly_qc = true;
    std::vector<int> block_weights(num_blocks_y * num_blocks_x, 0);
    
    for (int by = 0; by < num_blocks_y; by++) {
        for (int bx = 0; bx < num_blocks_x; bx++) {
            
            // Check first row of block
            int first_row_weight = 0;
            std::vector<int> first_row(block_w, 0);
            for (int c = 0; c < block_w; c++) {
                int r_global = by * block_h;
                int c_global = bx * block_w + c;
                if (img[r_global * width + c_global] < 128) {
                    first_row_weight++;
                    first_row[c] = 1;
                }
            }
            block_weights[by * num_blocks_x + bx] = first_row_weight;
            
            // Verify cyclic property
            for (int r = 1; r < block_h; r++) {
                for (int c = 0; c < block_w; c++) {
                    int r_global = by * block_h + r;
                    int c_global = bx * block_w + c;
                    int expected = first_row[(c - r + block_w) % block_w];
                    int actual = (img[r_global * width + c_global] < 128) ? 1 : 0;
                    if (actual != expected) {
                        is_perfectly_qc = false;
                    }
                }
            }
        }
    }
    
    std::cout << "Is perfectly Quasi-Cyclic? " << (is_perfectly_qc ? "True" : "False") << "\n";
    
    int global_min_block = block_w, global_max_block = 0;
    for (int w : block_weights) {
        if (w < global_min_block) global_min_block = w;
        if (w > global_max_block) global_max_block = w;
    }
    
    std::cout << "Block weights (first row of each block): min=" << global_min_block << ", max=" << global_max_block << "\n";
    
    // Print a map of block weights for the first 8 rows of blocks
    for (int by = 0; by < std::min(8, num_blocks_y); by++) {
        std::cout << "Block Row " << by << ": [";
        for (int bx = 0; bx < num_blocks_x; bx++) {
            std::cout << block_weights[by * num_blocks_x + bx] << (bx < num_blocks_x - 1 ? ", " : "");
        }
        std::cout << "]\n";
    }
    
    stbi_image_free(img);
    return 0;
}
