02-Core-Fitting-Pipeline

# LPGS Core Fitting Pipeline

> **Doc Version**: v1.0.0 · **Algorithm Version**: LPGS v1.0

## 1. Image Loading & Preprocessing

### 1.1 Resolution Normalization

The source image is first scaled so that its maximum dimension does not exceed `MAX_DIM = 500` pixels:

```
scale = min(1.0, MAX_DIM / max(srcWidth, srcHeight))
scaledW = round(srcWidth * scale)
scaledH = round(srcHeight * scale)
```

This step ensures that regardless of the input image size, fitting always operates within a fixed-order pixel space, guaranteeing predictable algorithm runtime.

### 1.2 Edge Padding (Edge-Aware Mode)

When edge-aware mode is enabled, a 3% transparent margin is added around the image:

```
edgePadding = ceil(max(scaledW, scaledH) * 0.03)
imgW = scaledW + edgePadding * 2
imgH = scaledH + edgePadding * 2
```

The padding provides "overflow space" for shapes near edges, preventing hard clipping. Target pixels within the padding area are set to fully transparent.

### 1.3 Alpha Pre-multiplication (Preserve Transparency)

When transparency preservation mode is enabled, the target image undergoes alpha pre-multiplication:

For each pixel:

- If `targetAlpha < alphaThreshold`: the pixel is treated as "transparent" — its target color is replaced with the current canvas color (equivalent to "ignoring this pixel")
- If `targetAlpha >= alphaThreshold`: semi-transparent pixels are composited with the background:
  ```
  targetRGB = targetRGB * (alpha/255) + currentRGB * (1 - alpha/255)
  targetAlpha = 255
  ```

After pre-multiplication, all pixels participating in fitting have Alpha = 255, converting the semi-transparent fitting problem into an opaque fitting problem.

### 1.4 Initial Error Calculation

After preprocessing, the initial total SSD is computed:

```
E₀ = Σ Σ (T_c(p) - C_c(p))²
     p  c∈{R,G,B}
```

where T is the target image, C is the current canvas (initialized to a solid background), and the summation covers only pixels above the Alpha threshold.

---

## 2. Error Grid

### 2.1 Grid Construction

The image is divided into `16×16` pixel grid cells, each recording the average SSD error of all pixels within it:

```
gridCell[gx, gy] = mean( SSD(p) for p in cell[gx, gy] )
```

The grid dimensions are `(imgW + 15) / 16 × (imgH + 15) / 16`, yielding approximately `31×31` cells for a 500px image.

### 2.2 Weighted Position Sampling

The error grid guides the placement of new shapes. The sampling process is as follows:

1. Compute the sum of all grid cell errors as `totalError`
2. Generate a random number `r` in the range `[0, totalError)`
3. Accumulate error values across grid cells in order; select the cell where the accumulated value first exceeds `r`
4. Uniformly sample a random coordinate point within the selected cell's pixel range

This is equivalent to **Roulette Wheel Selection** weighted by error values: regions with higher error have a higher probability of being selected.

### 2.3 Dirty Flag Mechanism

The error grid has a dirty flag `_errorGridDirty`. Whenever `ApplyShape()` modifies the current canvas, the flag is set to `true`, triggering a rebuild on the next sampling. This avoids the overhead of rebuilding the grid every frame.

---

## 3. Edge Distance Map

### 3.1 Construction Algorithm

The edge distance map is a float array of the same size as the image, recording the distance from each pixel to the nearest opaque edge. The construction uses a **Two-Pass Distance Transform**:

**Initialization**:
- Opaque pixels (`originalAlpha >= threshold`) → distance set to `maxDist = sqrt(imgW² + imgH²)`
- Transparent pixels → distance set to `0`

**Forward pass** (top-left → bottom-right):
```
For each pixel (x, y):
    minD = current value
    minD = min(minD, distMap[y-1, x]   + 1)      // up
    minD = min(minD, distMap[y, x-1]   + 1)      // left
    minD = min(minD, distMap[y-1, x-1] + 1.414)  // top-left
    minD = min(minD, distMap[y-1, x+1] + 1.414)  // top-right
    distMap[y, x] = minD
```

**Backward pass** (bottom-right → top-left):
```
For each pixel (x, y):
    minD = current value
    minD = min(minD, distMap[y+1, x]   + 1)      // down
    minD = min(minD, distMap[y, x+1]   + 1)      // right
    minD = min(minD, distMap[y+1, x+1] + 1.414)  // bottom-right
    minD = min(minD, distMap[y+1, x-1] + 1.414)  // bottom-left
    distMap[y, x] = minD
```

Both passes are executed twice (4 passes total) for a more accurate distance field.

### 3.2 Usage

The edge distance map serves two purposes in edge-aware mode:

1. **Limiting shape size**: The maximum radius of a shape is limited to 0.8 times the distance from its center to the nearest edge, preventing shapes from "overflowing" into transparent areas
2. **Overflow detection**: `ShapeOverlapsTransparent()` checks whether more than 5% of a shape's interior pixels fall within transparent areas

---

## 4. Candidate Shape Generation (RandomShape)

### 4.1 Position Sampling

- **Standard mode**: Position obtained via error grid weighted sampling
- **Precision mode**: Uniform random sampling within the ROI area

### 4.2 Shape Types & Parameters

LPGS supports three primitive types, each with a different parameter space:

| Type | Parameters | Position Sampling | Size Range | Rotation |
|------|-----------|-------------------|------------|----------|
| **circle** | X, Y, Radius | Error-weighted | `[4, 4 + minDim*0.2] × sizeMult` | None |
| **ellipse** | X, Y, Rx, Ry, Angle | Error-weighted | Rx: `[5, 5 + imgW*0.16] × sizeMult`<br>Ry: `[3, 3 + imgH*0.16] × sizeMult` | `[0, π)` uniform random |
| **line** | X, Y, Rx, Ry, Angle | Error-weighted | Rx: `[15, 15 + minDim*0.35] × sizeMult`<br>Ry: `[1.5, 4.5] × sizeMult` | Local gradient direction + π/2 |

### 4.3 Size Decay Strategy

Shape sizes decrease as fitting progresses:

```
progress = shapeCount / maxShapes
sizeMult = 1.0 - progress * 0.55
```

- Early fitting (`progress ≈ 0`): `sizeMult ≈ 1.0`, generating large shapes for the base
- Late fitting (`progress ≈ 1.0`): `sizeMult ≈ 0.45`, generating small shapes for detail

In precision mode, sizes are further reduced by 60%: `sizeMult *= 0.4`.

### 4.4 Line Gradient Alignment

The rotation angle of line primitives is not purely random, but based on the **local gradient direction** of the target image:

1. Take a 3×3 neighborhood around the sampling position
2. Compute the grayscale intensity of each pixel (BT.601 weights: `0.299R + 0.587G + 0.114B`)
3. Compute the gradient vector `(gx, gy)` using a Sobel-like approach
4. Line direction = gradient direction + π/2 (perpendicular to the gradient, i.e., along the edge direction)

If the gradient magnitude is too small (`gx² + gy² < 0.01`), it falls back to a random angle.

---

## 5. Shape Evaluation (EvaluateShape)

This is the mathematical core of LPGS — given a shape's geometric parameters, solve for its optimal color and opacity.

### 5.1 Problem Formulation

For the set of pixels S covered by the shape, the new pixel values after applying the shape are:

```
C'_c(p) = C_c(p) + α × (s_c - C_c(p)),  p ∈ S
```

where:
- `C_c(p)`: the value of color channel c at pixel p on the current canvas
- `s_c`: the shape's color channel c value (to be determined)
- `α`: the shape's opacity (to be determined)

The objective is to minimize the local SSD after applying the shape:

```
min_{s_R, s_G, s_B, α}  Σ_{p∈S} Σ_c (T_c(p) - C'_c(p))²
```

### 5.2 Two-Pass Iterative Least Squares

**First pass**:

1. Use the average color of target pixels as the initial estimate:

```
   s_c⁽⁰⁾ = (1/|S|) × Σ_{p∈S} T_c(p)
   ```

2. Given s_c⁽⁰⁾, solve for the optimal α:

```
   numerator = Σ_c [ s_c⁽⁰⁾ × (ΣT_c - ΣC_c) - (ΣT_cC_c - ΣC_c²) ]
   denominator = Σ_c [ |S| × (s_c⁽⁰⁾)² - 2 × s_c⁽⁰⁾ × ΣC_c + ΣC_c² ]
   α = numerator / denominator    (if denominator < 0.001, α defaults to 0.5)
   ```

3. Given α, solve for the optimal color:

```
   s_c = ( ΣT_c - (1-α) × ΣC_c ) / (α × |S|)
   ```

**Second pass**:

Recompute α using the s_c from the first pass, then recompute s_c. Two iterations are sufficient to converge to a high-precision solution.

### 5.3 Constraints

- `α ∈ [0.01, 0.99]`: avoid fully transparent or fully opaque shapes
- `s_c ∈ [0, 255]`: color values clamped to the valid range

### 5.4 Improvement Calculation

After evaluation, compute the shape's **improvement**:

```
oldSSD = covered area + uncovered area (before application)
newSSD = covered area (after application) + uncovered area (unchanged)
improvement = oldSSD - newSSD
newTotalSSD = totalError - oldSSD + newSSD
```

Note: `oldSSD` and `newSSD` are computed over the shape's bounding box, but only pixels inside the shape change — the error from outside pixels is identical in both and cancels out.

---

## 6. Candidate Selection & Mutation Optimization (TryOneShape)

### 6.1 Candidate Generation & Evaluation

1. Generate `Candidates` (default 64) random shapes
2. For each candidate, call `EvaluateShape()` to solve for optimal color
3. If edge-aware mode is enabled, check whether the shape overflows (`ShapeOverlapsTransparent`)
4. Select the candidate with the smallest `NewTotalSSD`

### 6.2 Mutation Optimization

Perform 16 rounds of random perturbation (mutation) on the best candidate. Each round:

1. Copy the current best shape
2. Apply random perturbation to position, size, and rotation:
   - Position: offset within `±0.3 × radius` range
   - Size: scale within `×(1 ± 0.3)` range
   - Rotation (ellipse): `±0.4 rad`, line: `±0.3 rad`
3. If edge-aware mode is enabled, constrain the mutated shape's size
4. Evaluate the mutated shape
5. If the mutation is better, replace the current best

### 6.3 Acceptance Criterion

A shape is accepted and applied to the canvas only when `bestSSD < totalError`. This guarantees monotonically decreasing total error.

---

## 7. Shape Application (ApplyShape)

Blend the shape's color into the current canvas:

```
For each pixel p covered by the shape:
    C_B(p) = round(C_B(p) + α × (s_B - C_B(p)))
    C_G(p) = round(C_G(p) + α × (s_G - C_G(p)))
    C_R(p) = round(C_R(p) + α × (s_R - C_R(p)))
```

Note that pixel storage format is BGRA (B comes first in the byte array), and the Alpha channel does not participate in blending.

After application, the error grid dirty flag is set to `true`, and the total error is updated to `bestSSD`.

---

## 8. Precision Fitting Mode

### 8.1 Reservation Mechanism

The standard fitting phase uses only `(1 - reserveRatio)` of the shape budget. For example, with `MaxShapes = 500, ReserveRatio = 20%`, standard fitting uses at most 400 shapes, leaving 100 for precision mode.

### 8.2 ROI Selection

Users can manually draw a rectangular region (ROI) on the output preview; precision fitting will only add shapes within that region.

### 8.3 Auto-Detection

`AutoDetectRoi()` slides a 1/4-canvas-sized window over the error grid and finds the region with the highest average error as the ROI.

### 8.4 Precision Mode Behavior

Upon entering precision mode:
- Shape positions are sampled only within the ROI
- Shape sizes are further reduced by 60% (`sizeMult *= 0.4`)
- `MaxShapes` is increased by `precisionAddNum` (default 50)