Tinkering with code and algorithms for fun
Public
Edited
Jun 11, 2023
2 stars
carved = {
const image = await FileAttachment("The_Scream_by_Edvard_Munch,_1893_-_Nasjonalgalleriet small@1.jpg").image();
const w = Math.min(width, image.naturalWidth);
const h = Math.round(w / image.naturalWidth * image.naturalHeight);
const context = DOM.context2d(w, h, 1);
context.canvas.style = "max-width: 100%;";
context.drawImage(image, 0, 0, w, h);
context.lineCap = "round";
context.lineJoin = "round";
context.lineWidth = 2;
const originalImageData = context.getImageData(0, 0, w, h);
const input = context.getImageData(0, 0, w, h);
const inputData = input.data;
const inputData32 = new Uint32Array(inputData.buffer);
// In order to put the removed seams in the right locations
// we have to store the original x-position of each pixel,
// so that we can re-insert them in the right order
const originalColumns = new Uint16Array(w * h);
const newColumns = new Uint16Array(w * h);
const removedSeamLoc = new Uint16Array(h);
const insertedSeamLoc = new Uint16Array(h);
function initBuffers() {
inputData.set(originalImageData.data);
for (let y = 0; y < h; y++) {
for (let x = 0; x < w; x++) {
originalColumns[x + y * w] = x;
}
}
newColumns.fill(0);
removedSeamLoc.fill(0);
insertedSeamLoc.fill(0);
}
const { sqrt, min, max, SQRT1_2 } = Math;
// ===============
// = Pixel delta =
// ===============
const pixelOrder = (function checkEndian() {
const arrayBuffer = new ArrayBuffer(2);
const uint8Array = new Uint8Array(arrayBuffer);
const uint16array = new Uint16Array(arrayBuffer);
uint8Array[0] = 0xAA; // set first byte
uint8Array[1] = 0xBB; // set second byte
if (uint16array[0] === 0xBBAA) return "ABGR"; // little endian
if (uint16array[0] === 0xAABB) return "RGBA"; // big endian
else throw new Error("Something crazy just happened");
})();
// This is probably slow compared to inlining it, but let's get it working first
const dPixel =
pixelOrder == "RGBA"
? function dPixel(p0, p1) {
const rgba0 = inputData32[p0];
const rgba1 = inputData32[p1];
const r0 = rgba0 >>> 24;
const g0 = (rgba0 >> 16) & 0xFF;
const b0 = (rgba0 >> 8) & 0xFF;
const r1 = rgba1 >>> 24;
const g1 = (rgba1 >> 16) & 0xFF;
const b1 = (rgba1 >> 8) & 0xFF;
// // Taking the difference between the squares
// // is a better approximation of linear rgb
// // while still being fairly cheap.
// const dr = r0*r0 - r1*r1;
// const dg = g0*g0 - g1*g1;
// const db = b0*b0 - b1*b1;
// return (dr < 0 ? -dr : dr) + (dg < 0 ? -dg : dg) + (db < 0 ? -db : db);
// Adapted from code in https://github.com/mapbox/pixelmatch
// calculate color difference according to the paper "Measuring perceived color difference
// using YIQ NTSC transmission color space in mobile applications" by Y. Kotsarenko and F. Ramos
const r = r0 - r1;
const g = g0 - g1;
const b = b0 - b1;
const y =
r * 0.2124681075446384 +
g * 0.4169973963260294 +
b * 0.08137907133969426;
const i =
r * 0.3258860837850668 -
g * 0.14992193838645426 -
b * 0.17596414539861255;
const q =
r * 0.0935501584120867 -
g * 0.23119531908149002 +
b * 0.13764516066940333;
return sqrt(y * y + i * i + q * q);
}
: function dPixel(p0, p1) {
const abgr0 = inputData32[p0];
const abgr1 = inputData32[p1];
const r0 = abgr0 & 0xFF;
const g0 = (abgr0 >> 8) & 0xFF;
const b0 = (abgr0 >> 16) & 0xFF;
const r1 = abgr1 & 0xFF;
const g1 = (abgr1 >> 8) & 0xFF;
const b1 = (abgr1 >> 16) & 0xFF;
// const dr = r0*r0 - r1*r1;
// const dg = g0*g0 - g1*g1;
// const db = b0*b0 - b1*b1;
// return (dr < 0 ? -dr : dr) + (dg < 0 ? -dg : dg) + (db < 0 ? -db : db);
const r = r0 - r1;
const g = g0 - g1;
const b = b0 - b1;
const y =
r * 0.2124681075446384 +
g * 0.4169973963260294 +
b * 0.08137907133969426;
const i =
r * 0.3258860837850668 -
g * 0.14992193838645426 -
b * 0.17596414539861255;
const q =
r * 0.0935501584120867 -
g * 0.23119531908149002 +
b * 0.13764516066940333;
return sqrt(y * y + i * i + q * q);
};
// =================
// = Cost Map Init =
// =================
// build up energy map
// see https://avikdas.com/2019/07/29/improved-seam-carving-with-forward-energy.html
const costMap = new Float64Array(w * h);
function initCostMap() {
// energy cost pixel top-left corner
costMap[0] = dPixel(0, 1);
// energy cost top row (except corners
for (let x = 1; x < w - 1; x++) {
costMap[x] = dPixel(x - 1, x + 1);
}
// energy cost pixel top-right corner
costMap[w - 1] = dPixel(w - 2, w - 1);
// build cost map for rest of the rows
for (let y = 1; y < h; y++) {
const line = y * w;
const linem1 = line - w;
// leftmost column
{
// const x = 0; // we can simply leave out the zero altogether here
const Cu = dPixel(line, line + 1) * 2;
const Cr = dPixel(linem1, line + 1);
const Mu = costMap[linem1] + Cu;
const Mr = costMap[linem1 + 1] + Cu + Cr;
costMap[line] = Mu < Mr ? Mu : Mr;
}
// rest of the row
for (let x = 1; x < w - 1; x++) {
const Cl = dPixel(x + linem1, x - 1 + line);
const Cu = dPixel(x - 1 + line, x + 1 + line);
const Cr = dPixel(x + linem1, x + 1 + line);
// Diagonals. We'll weight them sqrt(1/2) compared to the rest
const Clr = SQRT1_2 * dPixel(x - 1 + line, x + 1 + linem1); // bottom left to top right
const Crl = SQRT1_2 * dPixel(x + 1 + line, x - 1 + linem1); // bottom right to top left
const Ml = costMap[x - 1 + linem1] + Cu + Cl + Crl;
const Mu = costMap[x + linem1] + Cu + Clr + Crl;
const Mr = costMap[x + 1 + linem1] + Cu + Cr + Clr;
const M = Mr < Mu ? Mr : Mu;
costMap[x + line] = M < Ml ? M : Ml;
}
// rightmost column
{
const x = w - 1;
const Cl = dPixel(x + linem1, x - 1 + line);
const Cu = dPixel(x - 1 + line, x + line);
const Ml = costMap[x - 1 + linem1] + Cu + Cl;
const Mu = costMap[x + linem1] + Cu * 2;
costMap[x + line] = Ml < Mu ? Ml : Mu;
}
}
}
// ===============
// = Carve Loop =
// ===============
while (true) {
initBuffers();
initCostMap();
let inputEnd = w - 1;
while (inputEnd > 0) {
// find minimum of the bottom row
let y = h - 1;
let line = y * w;
let xMin = 0,
costMin = costMap[line];
for (let x = 1; x < inputEnd; x++) {
const xCost = costMap[x + line];
if (xCost < costMin) {
costMin = xCost;
xMin = x;
}
}
// remove seam - go from bottom minimum and carve up
while (y >= 0) {
// copy pixeldata + original location at xMin;
// moving 32 bits at once is faster than three 8-bit chunks
const rgbMin = inputData32[xMin + line];
const loc = originalColumns[xMin + line];
removedSeamLoc[y] = xMin;
// shift all remaining pixels (in "original" image) one to the left,
// including their locations
let i = xMin;
let idx = i + line;
while (i < inputEnd) {
inputData32[idx] = inputData32[idx + 1];
originalColumns[idx] = originalColumns[idx + 1];
idx++;
i++;
}
// insert removed pixel at the appropriate location on the right
while (i < w - 1 && loc > newColumns[idx + 1]) {
inputData32[idx] = inputData32[idx + 1];
newColumns[idx] = newColumns[idx + 1];
idx++;
i++;
}
// const idx4 = idx << 2;
inputData32[idx] = rgbMin;
newColumns[idx] = loc;
insertedSeamLoc[y] = i;
// go up by one row
y--;
line = y * w;
// The minimum cost is guaranteed to be one of the
// three pixels right above the pixel with the
// minimum cost of the row below it
costMin = costMap[xMin + line];
if (costMap[max(0, xMin - 1) + line] < costMin) {
xMin = max(0, xMin - 1);
costMin = costMap[xMin + line];
}
if (costMap[Math.min(inputEnd, xMin + 1) + line] < costMin) {
xMin = Math.min(inputEnd, xMin + 1);
costMin = costMap[xMin + line];
}
}
// Move the right border of the original (carved) image to the left by one
inputEnd--;
// Update image and draw seams
context.putImageData(input, 0, 0);
// context.strokeStyle = 'red';
// context.beginPath();
// context.moveTo(removedSeamLoc[0], 0);
// for (let y = 1; y < h; y++) {
// context.lineTo(removedSeamLoc[y], y);
// }
// context.stroke();
// context.strokeStyle = 'blue';
// context.beginPath();
// context.moveTo(insertedSeamLoc[0], 0);
// for (let y = 1; y < h; y++) {
// context.lineTo(insertedSeamLoc[y], y);
// }
// context.stroke();
yield {
context,
input,
inputData32,
inputEnd,
originalColumns,
newColumns
};
// Update cost map
// Optimization: we only have to update the pyramid
// under wherever xMin ends up being on the top row,
// because pixels outside of this pyramid cannot be
// affected by the change.
// Handle top-left and top-right corners
if (xMin === 0) {
costMap[0] = dPixel(0, 1);
costMap[1] = dPixel(0, 2);
} else if (xMin >= inputEnd) {
costMap[inputEnd - 1] = dPixel(inputEnd - 2, inputEnd);
costMap[inputEnd] = dPixel(inputEnd - 1, inputEnd);
} else {
costMap[xMin - 1] = dPixel(max(0, xMin - 2), xMin);
costMap[xMin] = dPixel(xMin - 1, xMin + 1);
costMap[xMin + 1] = dPixel(xMin, min(inputEnd, xMin + 2));
}
// update cost map for rest of the rows
let xStart = max(0, xMin - 1),
xEnd = min(xMin + 1, inputEnd);
for (let y = 1; y < h; y++) {
const line = y * w;
const linem1 = line - w;
let x = xStart;
// handle case where pyramid overlaps with start
if (x === 0) {
const Cu = dPixel(line, line + 1) * 2;
const Cr = dPixel(linem1, line + 1);
const Mu = costMap[linem1] + Cu;
const Mr = costMap[linem1 + 1] + Cu + Cr;
costMap[line] = Mu < Mr ? Mu : Mr;
} else {
const Cl = dPixel(x + linem1, x - 1 + line);
const Cu = dPixel(x - 1 + line, x + 1 + line);
const Cr = dPixel(x + linem1, x + 1 + line);
// Diagonals. We'll weight them sqrt(1/2) compared to the rest
const Clr = SQRT1_2 * dPixel(x - 1 + line, x + 1 + linem1); // bottom left to top right
const Crl = SQRT1_2 * dPixel(x + 1 + line, x - 1 + linem1); // bottom right to top left
const Ml = costMap[x - 1 + linem1] + Cu + Cl + Crl;
const Mu = costMap[x + linem1] + Cu + Clr + Crl;
const Mr = costMap[x + 1 + linem1] + Cu + Cr + Clr;
const M = Mr < Mu ? Mr : Mu;
costMap[x + line] = M < Ml ? M : Ml;
}
// inner pyramid
while (++x < xEnd) {
const Cl = dPixel(x + linem1, x - 1 + line);
const Cu = dPixel(x - 1 + line, x + 1 + line);
const Cr = dPixel(x + linem1, x + 1 + line);
const Ml = costMap[x - 1 + linem1] + Cu + Cl;
const Mu = costMap[x + linem1] + Cu;
const Mr = costMap[x + 1 + linem1] + Cu + Cr;
const M = Mr < Mu ? Mr : Mu;
costMap[x + line] = M < Ml ? M : Ml;
}
// handle case where pyramid overlaps with end
if (x === inputEnd) {
const Cl = dPixel(x + linem1, x - 1 + line);
const Cu = dPixel(x - 1 + line, x + line);
const Ml = costMap[x - 1 + linem1] + Cu + Cl;
const Mu = costMap[x + linem1] + Cu * 2;
costMap[x + line] = Ml < Mu ? Ml : Mu;
} else {
const Cl = dPixel(x + linem1, x - 1 + line);
const Cu = dPixel(x - 1 + line, x + 1 + line);
const Cr = dPixel(x + linem1, x + 1 + line);
const Ml = costMap[x - 1 + linem1] + Cu + Cl;
const Mu = costMap[x + linem1] + Cu;
const Mr = costMap[x + 1 + linem1] + Cu + Cr;
const M = Mr < Mu ? Mr : Mu;
costMap[x + line] = M < Ml ? M : Ml;
}
xStart = max(0, xStart - 1);
xEnd = min(inputEnd, xEnd + 1);
}
}
context.putImageData(input, 0, 0);
yield {
context,
input,
inputData32,
inputEnd,
originalColumns,
newColumns
};
}
}
One platform to build and deploy the best data apps
Experiment and prototype by building visualizations in live JavaScript notebooks. Collaborate with your team and decide which concepts to build out.
Use Observable Framework to build data apps locally. Use data loaders to build in any language or library, including Python, SQL, and R.
Seamlessly deploy to Observable. Test before you ship, use automatic deploy-on-commit, and ensure your projects are always up-to-date.