In an unrelated search, I just stumbled on this GitHub repo for optimal rectangle decomposition, using the "bipartite matching algorithm". It likely does a better job than what I've described below.
In the event that the link rots, here are the sources it cites for more information on the algorithm:
- Jr. W. Lipski, E. Lodi, F. Luccio, C. Mugnai, L. Pagli. (1979) "On two-dimensional data organization II". Fundamenta Informaticae
- E. Kawaguchi, T. Endo. (1980) "On a method of binary-picture representation and its application to data compression" IEEE Transactions on Pattern Analysis and Machine Intelligence
- T. Ohtsuki. (1982) "Minimum dissection of rectilinear regions" IEEE Conference on Circuits and Systems
- T. Suk, C. Hoschl, J. Flusser. (2012) "Rectangular Decomposition of Binary Images",Advanced Concepts for Intelligent Vision Systems. LNCS Vol. 7517
Older answer:
First, we can convert your source rectangles to cells in your underlying grid, to make the input more uniform. (Effectively rasterizing the problem)
This will let us find optimizations that might not be obvious when working directly with the source rectangles - particularly when it involves splitting multiple source rectangles to recombine them differently.
Next we can find connected regions of the same colour, using depth-first-search or flood filling algorithms. We can consider each connected region (a polyomino) in isolation - nothing that we do to a different region needs to influence this one.
Effectively we want to find a way to dissect this polyomino into rectangles (unfortunately most of the literature I can find is about the opposite problem: dissecting rectangles into polyominoes! This makes it tricky to search for leads...)
One straightforward method is to combine horizontal runs of adjacent squares into long skinny rectangles. Then we can compare against the row above and combine if our run starts & ends match up - either as we finish each run/row, or as we consider each cell to add to the current run.
I don't know yet how close this method gets to optimal. It seems it can run into a bit of trouble when a row it hasn't considered yet suggests a different split than the rows it's seen so far:
Detecting when a run/rectangle is exactly covered by runs above & below, then splitting it and merging them will solve this particular case, but I haven't explored how general the problem is.
I've also looked at methods where we walk the perimeter of the polyomino, and cut across anytime we encounter a concave corner, but this approach looks more error-prone to me. Getting optimal results seems to require prioritizing cuts that join two concave corners, and shapes containing hollows need special handling, so the row scan method seems to have the simplicity advantage.
One more method I'm looking at is to take the first run found in the top row & extend it down as far as you can. Then take the first run in the top row of what's left... This gets tripped-up on inverted T shapes though, so it's not optimal either.
I feel like there's probably a way to use dynamic programming to find the optimal split, but I haven't found it yet.