mHash — Marr–Hildreth Hash¶
Bring the image to 512×512, spread its gray levels evenly, and pass it through the Mexican-hat filter of Marr and Hildreth, which answers where the image has edges and fine detail at the filter's scale. Sum the answer over a 31×31 grid of blocks, and record, for each block of 64 small windows of that grid, whether it answers more strongly than its neighbors in the window. The result is 576 bits, a fingerprint of where the image's detail lies, much finer than the 64 bits of the other hashes and much more sensitive to where exactly it lies.
| Call | ph_compute_mhash(), or PH_ALGO_MHASH in ph_compute_digest() |
| Output | digest of 72 bytes, 576 bits, compared with ph_hamming_distance_digest() |
| Source | the operator: Marr and Hildreth, 1980; the hash built on it: pHash's ph_mh_imagehash(), which has no paper (provenance) |
| Cost | 1.8 ms on a 400×400 image, 91 ms on 20 Mpx, after decoding1 |
The steps¶
The steps below are those of the default parameters, a 512×512 image and a kernel reaching
8 pixels from its center;
ph_context_set_mhash_params() sets
others, below. Every picture on this page is computed by the library from
the same photograph, and the digest at the end is what ph_compute_mhash() returns for
it.
1. Grayscale. Each pixel becomes one luminance value,
the same integer approximation of the ITU-R BT.601 weights as for every grayscale hash.
ph_context_set_gray_weights() changes them; image preparation has why these weights.
2. Blur. The grayscale image is blurred at full resolution by a Gaussian of standard deviation 1 pixel, cut off at 3 pixels on either side and scaled so its seven weights add up to 1, first along the rows and then along the columns; the edges repeat their last pixel, and each result is rounded to a whole gray level.
3. Normalize to 512×512. The blurred image is resampled to 512×512 with the Mitchell filter, regardless of its aspect ratio: a wide picture is squeezed, a tall one stretched. Every image then reaches the filter at the same size, so the filter's scale is the same fraction of every picture.
4. Equalize. Each gray level \(v\) is replaced by its rank among the image's \(N = 512^2\) pixels,
where \(C(v)\) counts the pixels at level \(v\) or below and \(C_{\min}\) the pixels at the darkest level present. The darkest level becomes 0, the brightest 255, and the levels between spread so that each takes its share of the range.
The cumulative share becomes a straight line: whatever the exposure, the picture uses the whole range evenly. A change that keeps the order of the gray levels, such as a brighter exposure or a tone curve, changes the equalized image only through rounding.
5. The kernel. The filter is the Laplacian of Gaussian, sampled as the Mexican hat
with the defaults \(\alpha = 2\) and level \(\ell = 1\): \(s = 8\), a 17×17 kernel. Its center is positive, a ring around it negative, and its weights nearly cancel:
The 17×17 kernel at the defaults sums to 0.017, against a weight of 2 at its center.
The kernel crosses zero at \(A = 2\), a distance of \(\sqrt 2\, \alpha^{\ell} \approx 2.8\) pixels from its center. On a smooth area it answers almost nothing; on a line or a spot about that wide it answers strongly, positive on one side of a light-dark change and negative on the other.
6. The response. Each pixel of the equalized image \(I\) is replaced by the kernel's weighted sum around it, the edges repeating their last pixel:
The response is a map of detail: the spines of the cactus and the rim of the pot light up, the wall and the table stay pale.
7. Block sums. The response is summed over 16×16 blocks, 31 of them across and 31 down:
31 blocks of 16 cover 496 of the 512 pixels; the last 16 rows and columns enter only through the kernel's reach of the blocks next to them. The library does not compute \(R\) at every pixel; it computes these sums directly, faster and with the same result.
8. Windows and bits. The grid is read through 3×3 windows placed every 4 blocks, 8 across and 8 down: 64 windows, nine block sums each. The windows do not touch: every fourth row and column of the grid, seven of each, belongs to no window, so the digest reads 576 of the 961 block sums.
Each block sum in a window is compared with the mean of the nine, and sets its bit when it is above:
Each bit therefore says where the detail lies within its window, which sets the bits independently of the response's overall strength: the source scales the response to [0, 1] before this step, and the library leaves that out, since scaling every sum by the same positive factor and shifting it by the same amount cannot change which side of their mean any of them falls.
Folding the block sum into the kernel¶
Computing \(R\) at every pixel of the 496×496 area and then summing it costs \(496^2 \times 289\), some 71 million multiply-adds. The kernel weight does not depend on the pixel, so it can come out of the block sum:
and what is left inside is the sum of a 16×16 square of the image, shifted by \((j, i)\). An integral image, the table of sums of every top-left rectangle, gives any such square in four lookups. That makes \(31^2 \times 289\), some 280 thousand multiply-adds, 256 times fewer.
flowchart LR
I[Equalized image] --> R["Response at every pixel<br/>496² × 289"] --> S[Sum over each block] --> B[Block sums]
I --> T[Integral image] --> Q["289 shifted 16×16 squares<br/>per block, 31² × 289"] --> B
The square sums are exact integers, and only 289 products are added per block, in double precision. Evaluated the long way, the response sums 289 large products that nearly cancel, and in single precision the rounding error survives into the block sum. Most of the time it does not reach a bit; where it does, it is decisive. Here is the long way against the library's digest, with the response computed in double and in single precision:
| Images | Response computed in | Against the folded block sums |
|---|---|---|
| the example photograph | double | the same digest |
| the example photograph | float | the same digest |
| synthetic corpus, 24 images; not photographs | double | the same digest on every image |
| synthetic corpus, 24 images; not photographs | float | another digest on 1 of 24 images, 108 bits apart |
| 200 public-domain photographs from Wikimedia Commons | double | the same digest on every image |
| 200 public-domain photographs from Wikimedia Commons | float | another digest on 1 of 200 images, 1 bit apart |
The synthetic image where single precision decides 108 bits is four flat squares of color. In a window inside a flat area the nine block sums are equal, every comparison with their mean is a tie, and rounding alone breaks it. The folded sums break such ties the same way on every run; a direct evaluation in single precision breaks them by the order of its additions.
Digest layout¶
Bits run window by window, in the order of their numbers above, and within a window row by row: bit \(9w + 3r + c\) is row \(r\), column \(c\) of window \(w\). They are packed most significant bit first, eight to a byte, so the first window takes all of byte 0 and the top bit of byte 1.
The order is pHash's; neither paper specifies one. Two digests are compared bit by bit, so only the count of differing bits matters, not where they are.
What changes the hash¶
The hash records where fine detail lies, block by block, to within a few pixels. Anything that moves the detail across the grid changes it; anything that keeps the detail in place changes it little.
- Rotation and cropping move the hash furthest: a turn of a degree or two, or a crop of a few percent, moves the detail by a few pixels at the image's edges, as much as the kernel's center is wide, and the block sums near those edges change. A few degrees more, or a larger crop, and the digest is as far from the original as an unrelated image's, about half of its 576 bits.
- Brightness, contrast and gamma keep the order of the gray levels, which equalizing turns into the same image up to rounding, until the brightest pixels clip. They move a few percent of the bits.
- Recompression, noise, blur and downscaling change the fine detail itself, which is what the kernel measures; the hash moves with their strength, a few percent of the bits at moderate settings.
These are measurements on one photograph; the same edits over two whole corpora follow below.
The numbers behind the chart
| Transform | Strength | Bits that differ |
|---|---|---|
| JPEG quality | 95 | 1 |
| JPEG quality | 75 | 0 |
| JPEG quality | 50 | 35 |
| JPEG quality | 30 | 58 |
| JPEG quality | 15 | 69 |
| JPEG quality | 5 | 106 |
| Downscale | 0.75 | 17 |
| Downscale | 0.5 | 28 |
| Downscale | 0.25 | 90 |
| Downscale | 0.125 | 151 |
| Rotation | 1 | 106 |
| Rotation | 2 | 170 |
| Rotation | 5 | 268 |
| Rotation | 10 | 291 |
| Rotation | 20 | 280 |
| Rotation | 45 | 290 |
| Brightness | 0.5 | 26 |
| Brightness | 0.7 | 30 |
| Brightness | 0.85 | 23 |
| Brightness | 1.15 | 17 |
| Brightness | 1.3 | 56 |
| Brightness | 1.5 | 71 |
| Contrast | 0.5 | 26 |
| Contrast | 0.7 | 19 |
| Contrast | 0.85 | 17 |
| Contrast | 1.15 | 15 |
| Contrast | 1.3 | 24 |
| Contrast | 1.5 | 39 |
| Gamma | 0.5 | 30 |
| Gamma | 0.7 | 24 |
| Gamma | 0.85 | 22 |
| Gamma | 1.2 | 19 |
| Gamma | 1.5 | 17 |
| Gamma | 2 | 20 |
| Gaussian blur | 0.5 | 4 |
| Gaussian blur | 1 | 12 |
| Gaussian blur | 2 | 32 |
| Gaussian blur | 4 | 83 |
| Gaussian blur | 8 | 154 |
| Noise | 2 | 18 |
| Noise | 5 | 35 |
| Noise | 10 | 49 |
| Noise | 20 | 70 |
| Noise | 40 | 93 |
| Crop | 5 | 200 |
| Crop | 10 | 295 |
| Crop | 20 | 270 |
| Crop | 30 | 281 |
How this was measured
Each transform is applied alone to the photograph above, and the edited copy is
saved losslessly, as PPM (as a JPEG of the given quality for that panel). The library
hashes the original and every copy, and compares each copy with the original by
ph_similarity_digest(). The chart shows
\(576 \times (1 - \text{similarity})\), which is the Hamming distance: how many of the
576 bits the edit flipped. Below is the code that ran, not a copy of it.
The edits, in Python with Pillow:
def transforms():
"""(name, x-axis label, [(strength, PIL image -> (image, (file extension, quality)))]).
Each edit is applied alone to the original. A variant is saved as PPM, which is
lossless (and quick to write for a corpus of thousands of variants), so the edit is
the only difference; JPEG quality is the one edit that is the saving itself.
"""
def jpeg(q):
return lambda im: (im, ("jpg", q))
def scale(s):
return lambda im: (im.resize((max(1, round(im.width * s)), max(1, round(im.height * s))),
Image.LANCZOS), ("ppm", None))
def rotate(deg):
# The corners the turn uncovers are filled with the image's mean color, so they
# add no edge of their own.
def f(im):
fill = tuple(int(v) for v in np.asarray(im).reshape(-1, 3).mean(axis=0))
return im.rotate(deg, resample=Image.BICUBIC, fillcolor=fill), ("ppm", None)
return f
def brightness(k):
return lambda im: (ImageEnhance.Brightness(im).enhance(k), ("ppm", None))
def contrast(k):
return lambda im: (ImageEnhance.Contrast(im).enhance(k), ("ppm", None))
def gamma(g):
def f(im):
a = (np.asarray(im).astype(np.float64) / 255.0) ** g
return Image.fromarray((a * 255.0).round().astype(np.uint8)), ("ppm", None)
return f
def blur(r):
return lambda im: (im.filter(ImageFilter.GaussianBlur(r)), ("ppm", None))
def noise(sigma):
# A fixed seed: the same noise on every build, so the chart does not move.
def f(im):
rng = np.random.default_rng(1)
a = np.asarray(im).astype(np.float64) + rng.normal(0, sigma, (im.height, im.width, 3))
return Image.fromarray(np.clip(a, 0, 255).round().astype(np.uint8)), ("ppm", None)
return f
def crop(p):
def f(im):
dx, dy = round(im.width * p / 2), round(im.height * p / 2)
return im.crop((dx, dy, im.width - dx, im.height - dy)), ("ppm", None)
return f
return [
("JPEG quality", "quality", [(q, jpeg(q)) for q in (95, 75, 50, 30, 15, 5)]),
("Downscale", "scale", [(s, scale(s)) for s in (0.75, 0.5, 0.25, 0.125)]),
("Rotation", "degrees", [(d, rotate(d)) for d in (1, 2, 5, 10, 20, 45)]),
("Brightness", "factor", [(k, brightness(k)) for k in (0.5, 0.7, 0.85, 1.15, 1.3, 1.5)]),
("Contrast", "factor", [(k, contrast(k)) for k in (0.5, 0.7, 0.85, 1.15, 1.3, 1.5)]),
("Gamma", "exponent", [(g, gamma(g)) for g in (0.5, 0.7, 0.85, 1.2, 1.5, 2.0)]),
("Gaussian blur", "radius, px", [(r, blur(r)) for r in (0.5, 1, 2, 4, 8)]),
("Noise", "σ, gray levels", [(s, noise(s)) for s in (2, 5, 10, 20, 40)]),
("Crop", "border removed, %", [(p, crop(p / 100)) for p in (5, 10, 20, 30)]),
]
Writing the copies and handing them to the measuring tool:
def measure_variants(tool, base, ref):
"""{algorithm: {transform: [(strength, value)]}}: every edit of `base` (a Pillow
image), of both groups of transforms.py, compared with `ref`, the same pixels saved
losslessly.
Every variant is written to a file, and `site_stages measure` hashes them with the
library and compares each variant with the reference, by each algorithm's own
metric; `value` is None where the comparison does not apply.
"""
with tempfile.TemporaryDirectory() as tmp:
files, keys = [], []
for name, _, steps in transforms() + content_edits():
for strength, op in steps:
files.append(write_edit(op, base, os.path.join(tmp, str(len(files)))))
keys.append((name, strength))
rows = run_lines(tool, "measure", ref, *files)
result = {}
for (name, strength), row in zip(keys, rows):
for algo, value in row.items():
if algo != "file":
result.setdefault(algo, {}).setdefault(name, []).append((strength, value))
return result
def measure_robustness(tool, image):
"""measure_variants() for one image file."""
base = Image.open(image).convert("RGB")
with tempfile.TemporaryDirectory() as tmp:
ref = os.path.join(tmp, "reference.ppm")
base.save(ref)
return measure_variants(tool, base, ref)
The measuring tool, linked against the library: it computes every algorithm's digest of the original, then of each copy, and compares the two by that algorithm's own metric.
/* site_stages measure <reference> <variant>...: every algorithm's digest of the reference,
* then one JSON line per variant with its comparison to the reference, per algorithm (the
* similarity in [0, 1] for a bit hash, the peak correlation for Radial, the histogram
* intersection for ColorHash, the L2 distance for ColorMoments). A `--load=<settings>`
* between the images loads the ones after it so (take_load_settings()). */
int mode_measure(int argc, char **argv) {
ph_digest_t ref[PH_ALGORITHM_COUNT], d[PH_ALGORITHM_COUNT];
int have_ref[PH_ALGORITHM_COUNT], have[PH_ALGORITHM_COUNT];
int status = 0, images = 0;
for (int v = 0; v < argc && !status; v++) {
if (take_load_settings(argv[v], &status)) {
continue;
}
ph_context_t *ctx = NULL;
if (load_image(&ctx, argv[v])) {
return 1;
}
digest_all(ctx, images ? d : ref, images ? have : have_ref);
ph_free(ctx);
if (images++ == 0) {
continue;
}
json_t j = json_begin(stdout);
json_string(&j, "file", argv[v]);
compare_all(&j, ref, have_ref, d, have);
json_end(&j);
}
return status;
}
/* One algorithm's comparison of two digests by its own metric; 0 when it does not apply. */
static int compare(const ph_digest_t *a, const ph_digest_t *b, double *out) {
switch (a->kind) {
case PH_DIGEST_KIND_BITS:
*out = ph_similarity_digest(a, b);
return *out >= 0.0;
case PH_DIGEST_KIND_COEFFICIENTS:
return ph_radial_similarity(a, b, out) == PH_SUCCESS;
case PH_DIGEST_KIND_HISTOGRAM:
return ph_histogram_intersection(a, b, out) == PH_SUCCESS;
case PH_DIGEST_KIND_VECTOR16:
*out = ph_l2_distance(a, b);
return *out >= 0.0;
default:
return 0;
}
}
The number on the chart:
def bits_that_differ(similarity, bits_total):
"""The number on a bit hash's chart. ph_similarity_digest() is 1 - distance / bits,
so this is the Hamming distance itself: how many of the hash's bits the edit flipped."""
return (1.0 - similarity) * bits_total
To reproduce, make site builds the tool and redraws every figure. The tool can also
be run by hand on any pair of images; it prints one JSON line per variant, with each
algorithm's comparison:
Over two corpora¶
The same nine edits, applied to every image of two corpora: 200 public-domain photographs from Wikimedia Commons (the list), and the 24 generated images the library's tests measure (methodology), which are not photographs. Each line is the median over the corpus and the band around it holds the middle half of its images.
The photographs repeat the single image: rotation and cropping reach half of the bits at a few degrees and a few percent, the other edits stay within a small share of them. The synthetic images also react to noise. Most of them are flat areas of a few colors, whose narrow histograms equalizing stretches over the whole range, and with them any noise on the flat areas: the response then answers to the noise.
Copies and different images¶
Robustness alone proves little: a hash that never changes is perfectly robust. What makes a hash useful is that copies of an image land closer together than different images do. Here a copy is an original after one moderate edit (one strength of each of the nine, listed under "How this was measured"), and different images are every pair of distinct originals in the corpus.
Different images gather tightly around 288, half of the 576 bits: with that many bits, two unrelated digests differ in very nearly half of them. Copies split in two. Most land within a few dozen bits, but the rotated and cropped copies land near half, among the different images. The threshold that accepts 95 % of the copies has to reach up to them, and d′, the gap between the two means in units of their spread, is the lowest of the bit hashes on this site. That the threshold still lets few different pairs through is the tight cluster's doing, not the copies'.
The numbers behind the charts
| Corpus | Copies | Different pairs | d′ | Threshold | Different pairs within it |
|---|---|---|---|---|---|
| synthetic corpus, 24 images; not photographs | 216 | 276 | 2.94 | 255 | 6.9% |
| 200 public-domain photographs from Wikimedia Commons | 1800 | 19900 | 3.46 | 251 | 0.7% |
Median bits that differ over each corpus, with the 25th and 75th percentiles:
| Transform | Strength | synthetic corpus, 24 images; not photographs | 200 public-domain photographs from Wikimedia Commons |
|---|---|---|---|
| JPEG quality | 95 | 12.5 (2.5–48.5) | 2 (1–5) |
| JPEG quality | 75 | 39 (12.3–65) | 8 (4–20) |
| JPEG quality | 50 | 54.5 (15.2–91.5) | 12 (7–31.2) |
| JPEG quality | 30 | 70 (22–123.3) | 18 (9–42) |
| JPEG quality | 15 | 78 (1–169.5) | 31.5 (18.7–63) |
| JPEG quality | 5 | 108 (37.7–229.5) | 71.5 (47–108) |
| Downscale | 0.75 | 33 (20.5–54) | 7 (4–10.2) |
| Downscale | 0.5 | 68 (25.8–116) | 16 (11–23) |
| Downscale | 0.25 | 151 (59.5–224.5) | 48 (38.5–62) |
| Downscale | 0.125 | 208 (144.7–272.2) | 101 (86–120) |
| Rotation | 1 | 72.5 (30.5–124.3) | 144 (125–159.3) |
| Rotation | 2 | 111.5 (74.2–192) | 214 (191.5–228) |
| Rotation | 5 | 228.5 (124–258) | 273 (259–286) |
| Rotation | 10 | 257.5 (129.2–288.2) | 283 (272.7–293) |
| Rotation | 20 | 270.5 (186.7–287.3) | 286 (279–295) |
| Rotation | 45 | 278 (207.2–290.2) | 287 (276–293) |
| Brightness | 0.5 | 11 (0.7–40.7) | 7 (3–22) |
| Brightness | 0.7 | 8 (0–34.5) | 6 (2–16.3) |
| Brightness | 0.85 | 5.5 (0–31.5) | 5 (2–13.2) |
| Brightness | 1.15 | 6 (0–65) | 8 (3.7–20.2) |
| Brightness | 1.3 | 21 (0–134.7) | 14 (7–34) |
| Brightness | 1.5 | 35.5 (5.8–167.8) | 25.5 (13–48.5) |
| Contrast | 0.5 | 11 (0.7–49.3) | 7 (3–20.2) |
| Contrast | 0.7 | 10 (0–38.5) | 6 (2–17) |
| Contrast | 0.85 | 5.5 (0–32.5) | 5.5 (2–13.2) |
| Contrast | 1.15 | 7.5 (0–23.5) | 8 (3–18) |
| Contrast | 1.3 | 12 (0–57) | 12 (5.8–24.3) |
| Contrast | 1.5 | 27 (0.7–115.5) | 18.5 (8–31) |
| Gamma | 0.5 | 13 (0–135.8) | 15 (9–25) |
| Gamma | 0.7 | 6.5 (0–97.5) | 10 (5–18) |
| Gamma | 0.85 | 3.5 (0–69.3) | 6 (3–16) |
| Gamma | 1.2 | 3.5 (0–75.5) | 7 (4–14.3) |
| Gamma | 1.5 | 9 (0–153.2) | 12.5 (8–20) |
| Gamma | 2 | 22 (0–180.2) | 19 (12.7–30.2) |
| Gaussian blur | 0.5 | 4.5 (0–15.8) | 3 (1–5) |
| Gaussian blur | 1 | 25.5 (10–43.8) | 7 (4–12) |
| Gaussian blur | 2 | 83.5 (29.2–125) | 18.5 (13–24.3) |
| Gaussian blur | 4 | 169 (65.3–242) | 46 (37.5–59) |
| Gaussian blur | 8 | 208.5 (152.2–273.5) | 99 (82–114.3) |
| Noise | 2 | 91.5 (17.2–205.5) | 7 (3–24.3) |
| Noise | 5 | 132.5 (34.5–235.5) | 11 (5–36) |
| Noise | 10 | 170 (60.5–254) | 17 (9–52.3) |
| Noise | 20 | 196.5 (109–259.3) | 31 (17.8–74.7) |
| Noise | 40 | 237 (148.5–265.5) | 53 (34–102) |
| Crop | 5 | 173.5 (94.3–230.3) | 247 (229.3–259) |
| Crop | 10 | 240.5 (136.5–277.5) | 286 (269.5–298) |
| Crop | 20 | 275 (157.3–305.3) | 284.5 (273–293) |
| Crop | 30 | 265 (178.2–288) | 286 (276–295) |
How this was measured
For each corpus, every image goes through the nine edits above, and each copy is
compared with its original exactly as on the single photograph. Every pair of
distinct originals is compared too, by site_stages pairs. Below is the code that
ran.
The photographs are downloaded once and checked against the SHA-256 the manifest records for each:
def download(row, unreachable):
"""The photograph's bytes as the manifest names them, from the copy of the corpus or
else from Commons; (None, why) when neither has them. A source that cannot be reached
is added to `unreachable` and not tried again."""
why = "no source"
for source in (MIRROR + row["sha256"] + extension(row), row["url"]):
host = urllib.parse.urlsplit(source).netloc
if host in unreachable:
continue
try:
data = get(source)
except urllib.error.HTTPError as e:
why = f"{host}: {e}"
continue
except OSError as e:
unreachable[host] = str(e)
why = f"cannot reach {host} ({e})"
continue
finally:
if source == row["url"]:
time.sleep(0.2) # Commons asks clients to pace their requests
if hashlib.sha256(data).hexdigest() == row["sha256"]:
return data, None
why = f"{host}: the download differs from the manifest"
return None, why
def unpack_archive(rows, d):
"""Writes the rows' photographs that the copy's archive holds, verified, into d.
An archive that cannot be downloaded or read writes nothing; the files are then
downloaded one by one."""
wanted = {row["sha256"] + extension(row): row for row in rows}
try:
req = urllib.request.Request(MIRROR_ARCHIVE, headers={"User-Agent": USER_AGENT})
with urllib.request.urlopen(req, timeout=300) as r, \
tarfile.open(fileobj=r, mode="r|gz") as tar:
for member in tar:
row = wanted.get(os.path.basename(member.name))
if row is None or not member.isfile():
continue
data = tar.extractfile(member).read()
if hashlib.sha256(data).hexdigest() != row["sha256"]:
continue
path = os.path.join(d, local_name(row))
with open(path + ".part", "wb") as f:
f.write(data)
os.replace(path + ".part", path)
except (OSError, tarfile.TarError) as e:
print(f"warning: photo corpus: the archive of the copy: {e}", file=sys.stderr)
_fetched = {}
def fetch(cache=None):
"""The verified local paths of the manifest's photographs, downloading what is missing.
Returns (paths, problems); a photograph that is unavailable is a problem, not an error.
One process downloads once: a second call returns the first one's answer."""
d = cache_dir(cache)
if d in _fetched:
return _fetched[d]
os.makedirs(d, exist_ok=True)
def cached(row):
path = os.path.join(d, local_name(row))
if not os.path.exists(path):
return False
with open(path, "rb") as f:
return hashlib.sha256(f.read()).hexdigest() == row["sha256"]
rows = read_manifest()
missing = [row for row in rows if not cached(row)]
if len(missing) > 10:
unpack_archive(missing, d)
paths, problems, unreachable = [], [], {}
for row in rows:
path = os.path.join(d, local_name(row))
if row not in missing or cached(row):
paths.append(path)
continue
data, why = download(row, unreachable)
if data is None:
problems.append(f"{row['file']}: {why}")
continue
with open(path + ".part", "wb") as f:
f.write(data)
os.replace(path + ".part", path)
paths.append(path)
_fetched[d] = (paths, problems)
return paths, problems
Measuring a corpus, or taking the result from the cache when nothing that decides it has changed:
def measure_corpus(tool, corpus):
"""{"n", "label", "robust": {algo: {transform: [[strength, [value per image]]]}},
"edits": {the same, for the edits that change the picture}, "different": {algo:
[value per pair]}} for one corpus, from the cache when nothing that decides it has
changed. Values are each algorithm's own metric, as `site_stages measure` and
`site_stages pairs` print it; None where it does not apply."""
paths = images(tool, corpus)
def measure():
with tempfile.TemporaryDirectory() as tmp:
# The originals are decoded once, by Pillow, and written losslessly, so that
# the reference an edit is compared with and the image that is edited are the
# same pixels.
refs = [os.path.join(tmp, f"{i:04d}.ppm") for i in range(len(paths))]
with concurrent.futures.ProcessPoolExecutor() as pool:
per_image = list(pool.map(_measure_one, [(tool, p, r) for p, r in
zip(paths, refs)]))
rows = run_lines(tool, "pairs", *refs) if len(refs) >= 2 else []
robust, edits = {}, {}
content = {name for name, _, _ in content_edits()}
for result in per_image:
for algo, by_transform in result.items():
for name, points in by_transform.items():
group = edits if name in content else robust
slot = group.setdefault(algo, {}).setdefault(
name, [[s, []] for s, _ in points])
for (strength, value), cell in zip(points, slot):
cell[1].append(value)
different = {}
for row in rows:
for algo, value in row.items():
if algo not in ("a", "b"):
different.setdefault(algo, []).append(value)
return {"n": len(paths), "label": label(corpus, len(paths)), "robust": robust,
"edits": edits, "different": different}
return cached(corpus, cache_key(paths), measure, f"{corpus} corpus, {len(paths)} images")
Comparing every pair of originals:
/* site_stages pairs <image>...: every image's digests, then one JSON line per pair of
* distinct images (i < j, in the order given, as "a" and "b") with their comparison,
* per algorithm. A `--load=<settings>` before the images loads them so. */
int mode_pairs(int argc, char **argv) {
ph_digest_t(*dig)[PH_ALGORITHM_COUNT] = calloc((size_t)argc, sizeof(*dig));
int (*have)[PH_ALGORITHM_COUNT] = calloc((size_t)argc, sizeof(*have));
int bad = (!dig || !have) && fail("out of memory", NULL);
int n = 0;
for (int i = 0; i < argc && !bad; i++) {
if (take_load_settings(argv[i], &bad)) {
continue;
}
ph_context_t *ctx = NULL;
bad = load_image(&ctx, argv[i]);
if (!bad) {
digest_all(ctx, dig[n], have[n]);
n++;
ph_free(ctx);
}
}
for (int i = 0; i < n && !bad; i++) {
for (int k = i + 1; k < n; k++) {
json_t j = json_begin(stdout);
json_int(&j, "a", i);
json_int(&j, "b", k);
compare_all(&j, dig[i], have[i], dig[k], have[k]);
json_end(&j);
}
}
free(dig);
free(have);
return bad;
}
The edits that make a copy, and the numbers on the second chart:
# A "copy" is an original after one moderate edit: one strength of each transform, the
# kind of change a picture goes through between two places it is stored.
COPY_STRENGTHS = {
"JPEG quality": 50,
"Downscale": 0.5,
"Rotation": 2,
"Brightness": 1.15,
"Contrast": 1.15,
"Gamma": 1.2,
"Gaussian blur": 1,
"Noise": 5,
"Crop": 5,
}
def separability(copies, different, lower_is_closer, recall=RECALL):
"""d', the threshold that keeps `recall` of the copies (RECALL unless a page asks for
another share), and the share of different pairs that threshold also accepts.
d' is the gap between the two means in units of their pooled spread, as
tests/src/test_hash_properties.c defines it (separability()), signed so that a
positive d' means copies are closer than different images.
"""
c, d = np.asarray(copies, float), np.asarray(different, float)
sign = 1.0 if lower_is_closer else -1.0
pooled = np.sqrt((c.std() ** 2 + d.std() ** 2) / 2.0)
dprime = sign * (d.mean() - c.mean()) / pooled if pooled > 0 else float("inf")
# The threshold: the closest value that still accepts `recall` of the copies.
t = sign * np.quantile(sign * c, recall, method="inverted_cdf")
accepted = (d <= t) if lower_is_closer else (d >= t)
return dprime, t, accepted.mean()
The synthetic corpus is generated by make_base() in
tests/src/synthetic_corpus.h, the same code the tests use. make site measures
both corpora and draws every figure.
Edits that change the picture¶
The edits above are those a copy goes through. The ones here change what the picture shows: a patch of the image's own top-left corner pasted over its center, covering 1 to 16 % of the frame; the hue of every pixel turned, which recolors the picture; and a quarter turn, a half turn and a mirror image. Whether a hash should notice them depends on what it is for: a search for copies wants to ignore a recoloring, a search for retouched pictures wants to catch the patch.
Each panel shows how far the edited images land from their originals over both corpora. The dashed line is the threshold of the chart above, the one that accepts 95 % of the copies, in each corpus's color: an edited image below it would be taken for a copy.
- Every patch passes for a copy. The threshold sits so high, to take in the rotated and cropped copies, that even a patch over 16 % of the frame lands far below it.
- Turning the hue moves the hash little on photographs, whose muted colors keep their luminance under a hue rotation, and far on the flat, saturated colors of the synthetic images. Equalizing does not protect against it: a recoloring changes the order of the gray levels, not only their spacing.
- A turn or a mirror image moves the hash as far as an unrelated image.
The numbers behind the chart
Median bits that differ, and the share of edited images within the threshold that accepts 95% of the copies; synthetic: synthetic corpus, 24 images; not photographs; photographs: 200 public-domain photographs from Wikimedia Commons.
| Edit | Synthetic: median (25th–75th) | accepted as a copy | Photographs: median (25th–75th) | accepted as a copy |
|---|---|---|---|---|
| patch over 1% of the frame | 9.5 (2.5–22.2) | 100% | 5 (3–8) | 100% |
| patch over 4% of the frame | 21.5 (17.5–29.2) | 100% | 19 (16.8–23) | 100% |
| patch over 9% of the frame | 47 (20.7–60.8) | 100% | 37 (33–44) | 100% |
| patch over 16% of the frame | 59.5 (49.3–68.5) | 100% | 59 (53.8–64) | 100% |
| hue turned 30° | 67 (4–225.5) | 88% | 31 (18–48) | 100% |
| hue turned 90° | 268.5 (87.7–429.8) | 46% | 25 (13–37) | 100% |
| hue turned 180° | 298 (106–402.8) | 42% | 49 (28.7–75.2) | 100% |
| turned 90° | 287.5 (119.5–297.2) | 38% | 284 (273–293) | 8% |
| turned 180° | 218.5 (82.2–328.5) | 54% | 288 (274.8–298.2) | 10% |
| mirrored | 273.5 (125–360.2) | 46% | 278 (257–294) | 22% |
How this was measured
The same measurement as for the copies, over the same images: each edit is applied alone to the original, the copy is saved losslessly, and the library compares it with the original. The threshold is the one of the chart above, computed from the copies. Below is the code of the edits.
def content_edits():
"""The edits that change the picture, in the shape of transforms(); a strength is a
number, or a word where the edit has no scale.
A patch replaces a square at the center, covering the given share of the frame, with
the square of the same size at the top-left corner of the same image: a local edit
with content the image itself supplies. A hue rotation turns every pixel's hue by the
given angle and keeps its saturation and value, a recoloring. The turns and the mirror
change only the orientation.
"""
def patch(area):
def f(im):
side = area ** 0.5
w, h = round(im.width * side), round(im.height * side)
x0, y0 = (im.width - w) // 2, (im.height - h) // 2
out = im.copy()
out.paste(im.crop((0, 0, w, h)), (x0, y0))
return out, ("ppm", None)
return f
def hue(degrees):
def f(im):
hsv = np.asarray(im.convert("HSV")).copy()
hsv[..., 0] = (hsv[..., 0].astype(int) + round(degrees * 256 / 360)) % 256
return Image.fromarray(hsv, "HSV").convert("RGB"), ("ppm", None)
return f
def turn(method):
return lambda im: (im.transpose(method), ("ppm", None))
return [
("Local patch", "share of the frame, %", [(a, patch(a / 100)) for a in (1, 4, 9, 16)]),
("Hue rotation", "degrees", [(d, hue(d)) for d in (30, 90, 180)]),
("Turn and mirror", "", [("90°", turn(Image.Transpose.ROTATE_90)),
("180°", turn(Image.Transpose.ROTATE_180)),
("mirror", turn(Image.Transpose.FLIP_LEFT_RIGHT))]),
]
A small edit against a rescale¶
A hash that should find edited copies of a picture has to move further for an edit than for the changes a copy goes through. Here are mHash's median distances for both, the edits of a copy at the strength that makes one, and the patches:
On the photographs, a patch over 4 % of the frame moves mHash about as far as halving the image's size, and a 2-degree turn or a 5 % crop moves it several times further than a patch over 16 %. On the synthetic images, recompression, noise and downscaling move it further than small patches. The order is the wrong way round for finding edits: mHash describes where detail lies, and any small shift of the detail outweighs a local change of it. For aHash, dHash, pHash and wHash the order is the right one: halving the size moves none of their bits, a patch moves some.
How the nine algorithms compare is on choosing an algorithm.
Cost¶
mHash is the most expensive of the hashes on a small image. The grayscale conversion, the blur and the resampling read every pixel of the image and grow with it; the equalizing and the integral image work at 512×512 whatever the image; the block sums take a fixed 280 thousand multiply-adds. On a large image the blur, at full resolution, is most of the time.
How this was measured
Measured on Intel(R) Xeon(R) Platinum 8370C CPU @ 2.80GHz, 4 cores, Linux x86_64; compiled by 13.3.0, optimized, with libjpeg-turbo decoding the JPEGs; the minimum of the runs. Decoding, every hash at its default settings, and in bold this page's.
| Case | 400×400 | 20 Mpx |
|---|---|---|
decode |
0.53 ms | 105 ms |
ahash |
0.21 ms | 23 ms |
dhash |
0.23 ms | 27 ms |
phash |
0.21 ms | 23 ms |
whash |
0.21 ms | 23 ms |
bmh |
0.21 ms | 23 ms |
mhash |
1.8 ms | 91 ms |
radial |
1.4 ms | 120 ms |
color_hash |
0.24 ms | 31 ms |
color_moments |
0.52 ms | 66 ms |
mhash_size_62 |
0.85 ms | 92 ms |
mhash_size_128 |
0.91 ms | 94 ms |
mhash_size_256 |
1.1 ms | 93 ms |
mhash_size_512 |
1.8 ms | 92 ms |
mhash_size_1024 |
4.7 ms | 97 ms |
mhash_size_2048 |
17 ms | 109 ms |
mhash_size_4096 |
66 ms | 159 ms |
Each case runs once to warm up, then until it has run at least five times and for at least a second; the time is the minimum. The grayscale image and the area grid the library caches are dropped before every run, so each run is the first hash on a loaded image. The cases, and the timing loop:
/* Seconds on a monotonic clock, as tests/src/bench_hash.c reads it. The calendar clock
* (timespec_get() with TIME_UTC) would step back when the system corrects its time, and a
* run spanning the step would come out negative and become the minimum. On Darwin,
* CLOCK_MONOTONIC counts whole microseconds, 2 % of the cheapest hash; mach_absolute_time()
* counts nanoseconds. */
static double seconds(void) {
#ifdef __APPLE__
static mach_timebase_info_data_t tb;
if (tb.denom == 0) {
mach_timebase_info(&tb);
}
return (double)mach_absolute_time() * tb.numer / tb.denom * 1e-9;
#else
struct timespec ts;
clock_gettime(CLOCK_MONOTONIC, &ts);
return (double)ts.tv_sec + (double)ts.tv_nsec * 1e-9;
#endif
}
/* One timed case: `run` does the work once on a loaded image and returns nonzero on
* failure; `arg` selects a variant of it. A case named after an algorithm, or starting
* with its name and "_", is shown on that algorithm's page. */
typedef struct {
const char *name;
int (*run)(ph_context_t *ctx, const char *path, int arg);
int arg;
} time_case_t;
static int time_decode(ph_context_t *ctx, const char *path, int arg) {
(void)arg;
return ph_load_from_file(ctx, path) != PH_SUCCESS;
}
/* The first hash computed on a loaded image, which is what a caller who loads an image
* and hashes it pays: the grayscale conversion and the shared area grid are dropped
* before every run, so each run computes them (the dropping is untimed). */
static int time_hash(ph_context_t *ctx, const char *path, int algo) {
(void)path;
ph_digest_t d;
return ph_compute_digest(ctx, (ph_algorithm_t)algo, &d) != PH_SUCCESS;
}
/* The four 64-bit hashes in one ph_compute_multi() call, which computes the grayscale
* image and the area grid once for all of them. */
static int time_multi(ph_context_t *ctx, const char *path, int arg) {
(void)path;
(void)arg;
uint64_t h[PH_HASH_FLAGS_COUNT];
return ph_compute_multi(ctx, PH_HASH_AHASH | PH_HASH_DHASH | PH_HASH_PHASH | PH_HASH_WHASH,
h) != PH_SUCCESS;
}
static int time_whash_full(ph_context_t *ctx, const char *path, int arg) {
(void)path;
(void)arg;
uint64_t h;
int bad = ph_context_set_whash_mode(ctx, PH_WHASH_FULL) != PH_SUCCESS ||
ph_compute_whash(ctx, &h) != PH_SUCCESS;
return ph_context_set_whash_mode(ctx, PH_WHASH_FAST) != PH_SUCCESS || bad;
}
static int time_mhash_size(ph_context_t *ctx, const char *path, int size) {
(void)path;
ph_digest_t d;
int bad = ph_context_set_mhash_params(ctx, PH_MH_ALPHA, PH_MH_LEVEL, size) != PH_SUCCESS ||
ph_compute_mhash(ctx, &d) != PH_SUCCESS;
return ph_context_set_mhash_params(ctx, PH_MH_ALPHA, PH_MH_LEVEL, PH_MH_IMAGE_SIZE) !=
PH_SUCCESS ||
bad;
}
/* pHash with a dct_size other than the default, and the default block. */
static int time_phash_dct(ph_context_t *ctx, const char *path, int size) {
(void)path;
uint64_t h;
int bad = ph_context_set_phash_params(ctx, size, PH_DCT_REDUCTION_SIZE) != PH_SUCCESS ||
ph_compute_phash(ctx, &h) != PH_SUCCESS;
return ph_context_set_phash_params(ctx, PH_DCT_SIZE, PH_DCT_REDUCTION_SIZE) != PH_SUCCESS ||
bad;
}
/* Radial with the settings of radial_variants[variant]. */
static int time_radial(ph_context_t *ctx, const char *path, int variant) {
(void)path;
ph_digest_t d;
int bad = !radial_settings(ctx, &radial_variants[variant]) ||
ph_compute_radial_hash(ctx, &d) != PH_SUCCESS;
return !radial_settings(ctx, &radial_variants[RADIAL_DEFAULT]) || bad;
}
/* Loading the image, then, unless `algo` is LOAD_ONLY, computing that algorithm: what a
* caller pays for a load and a hash. `arg` packs the load -- a ph_decode_scale_t, plus
* LOAD_GRAY for ph_context_set_load_grayscale() and LOAD_FRESH for a context created and
* freed around the load, as a caller who creates one per image does -- and the algorithm
* above it. Without LOAD_FRESH, a context of its own per kind of load, kept across runs,
* so the cases after it hash the image as it was first loaded. */
enum {
LOAD_GRAY = 4,
LOAD_FRESH = 8,
LOAD_ALGO = 16,
LOAD_ONLY = -1,
};
static int time_load(ph_context_t *ctx, const char *path, int arg) {
(void)ctx;
static ph_context_t *kept[LOAD_FRESH];
int kind = arg % LOAD_FRESH, fresh = arg & LOAD_FRESH, algo = arg / LOAD_ALGO - 1;
ph_context_t *c = fresh ? NULL : kept[kind];
if (!c && ph_create(&c) != PH_SUCCESS) {
return 1;
}
ph_digest_t d;
int bad = ph_context_set_decode_scale(c, (ph_decode_scale_t)(kind & 3)) != PH_SUCCESS ||
ph_context_set_load_grayscale(c, (kind & LOAD_GRAY) != 0) != PH_SUCCESS ||
ph_load_from_file(c, path) != PH_SUCCESS ||
(algo != LOAD_ONLY && ph_compute_digest(c, (ph_algorithm_t)algo, &d) != PH_SUCCESS);
if (fresh) {
ph_free(c);
} else {
kept[kind] = c;
}
return bad;
}
#define LOADED(name, flags, algo) {name, time_load, (flags) + ((algo) + 1) * LOAD_ALGO}
static const time_case_t time_cases[] = {
{"decode", time_decode, 0},
{"ahash", time_hash, PH_ALGO_AHASH},
{"dhash", time_hash, PH_ALGO_DHASH},
{"phash", time_hash, PH_ALGO_PHASH},
{"whash", time_hash, PH_ALGO_WHASH},
{"whash_full", time_whash_full, 0},
{"bmh", time_hash, PH_ALGO_BMH},
{"mhash", time_hash, PH_ALGO_MHASH},
{"radial", time_hash, PH_ALGO_RADIAL},
{"color_hash", time_hash, PH_ALGO_COLOR_HASH},
{"color_moments", time_hash, PH_ALGO_COLOR_MOMENTS},
{"multi", time_multi, 0},
{"mhash_size_62", time_mhash_size, 62},
{"mhash_size_128", time_mhash_size, 128},
{"mhash_size_256", time_mhash_size, 256},
{"mhash_size_512", time_mhash_size, 512},
{"mhash_size_1024", time_mhash_size, 1024},
{"mhash_size_2048", time_mhash_size, 2048},
{"mhash_size_4096", time_mhash_size, 4096},
{"phash_dct_8", time_phash_dct, 8},
{"phash_dct_16", time_phash_dct, 16},
{"phash_dct_24", time_phash_dct, 24},
{"radial_sigma_1", time_radial, RADIAL_SIGMA_1},
{"radial_sigma_8", time_radial, RADIAL_SIGMA_8},
{"radial_gamma_2", time_radial, RADIAL_GAMMA_2},
{"radial_grid_40x32", time_radial, RADIAL_GRID_40X32},
{"radial_grid_90x64", time_radial, RADIAL_GRID_90X64},
{"radial_grid_360x256", time_radial, RADIAL_GRID_360X256},
{"radial_grid_1440x1024", time_radial, RADIAL_GRID_1440X1024},
{"radial_grid_4096x4096", time_radial, RADIAL_GRID_4096X4096},
LOADED("scale_half_decode", PH_DECODE_SCALE_HALF, LOAD_ONLY),
LOADED("scale_quarter_decode", PH_DECODE_SCALE_QUARTER, LOAD_ONLY),
LOADED("scale_eighth_decode", PH_DECODE_SCALE_EIGHTH, LOAD_ONLY),
LOADED("scale_full_phash", PH_DECODE_SCALE_FULL, PH_ALGO_PHASH),
LOADED("scale_half_phash", PH_DECODE_SCALE_HALF, PH_ALGO_PHASH),
LOADED("scale_quarter_phash", PH_DECODE_SCALE_QUARTER, PH_ALGO_PHASH),
LOADED("scale_eighth_phash", PH_DECODE_SCALE_EIGHTH, PH_ALGO_PHASH),
LOADED("scale_full_mhash", PH_DECODE_SCALE_FULL, PH_ALGO_MHASH),
LOADED("scale_half_mhash", PH_DECODE_SCALE_HALF, PH_ALGO_MHASH),
LOADED("scale_quarter_mhash", PH_DECODE_SCALE_QUARTER, PH_ALGO_MHASH),
LOADED("scale_eighth_mhash", PH_DECODE_SCALE_EIGHTH, PH_ALGO_MHASH),
LOADED("scale_full_radial", PH_DECODE_SCALE_FULL, PH_ALGO_RADIAL),
LOADED("scale_half_radial", PH_DECODE_SCALE_HALF, PH_ALGO_RADIAL),
LOADED("scale_quarter_radial", PH_DECODE_SCALE_QUARTER, PH_ALGO_RADIAL),
LOADED("scale_eighth_radial", PH_DECODE_SCALE_EIGHTH, PH_ALGO_RADIAL),
LOADED("gray_decode", LOAD_GRAY, LOAD_ONLY),
LOADED("gray_phash", LOAD_GRAY, PH_ALGO_PHASH),
LOADED("gray_mhash", LOAD_GRAY, PH_ALGO_MHASH),
LOADED("gray_radial", LOAD_GRAY, PH_ALGO_RADIAL),
LOADED("fresh_decode", LOAD_FRESH, LOAD_ONLY),
LOADED("fresh_phash", LOAD_FRESH, PH_ALGO_PHASH),
};
static int cmp_double(const void *a, const void *b) {
double x = *(const double *)a, y = *(const double *)b;
return (x > y) - (x < y);
}
/* One warm-up run, then at least TIME_MIN_RUNS runs and TIME_MIN_SECONDS, at most
* TIME_MAX_RUNS; the minimum is what the pages quote. */
enum {
TIME_MIN_RUNS = 5,
TIME_MAX_RUNS = 300,
};
#define TIME_MIN_SECONDS 1.0
/* Whether case `name` is among the `n` names of `only`; every case is when `n` is 0. */
static int time_chosen(const char *name, int n, char **only) {
for (int i = 0; i < n; i++) {
if (strcmp(name, only[i]) == 0) {
return 1;
}
}
return n == 0;
}
/* site_stages time <image> [case...]: every case of time_cases[] on the image, or the cases
* named; the minimum and the median of its runs, in milliseconds, as one JSON object, with
* the build that ran them. */
int mode_time(int argc, char **argv) {
for (int i = 1; i < argc; i++) {
int known = 0;
for (size_t c = 0; c < COUNT(time_cases); c++) {
known |= strcmp(argv[i], time_cases[c].name) == 0;
}
if (!known) {
return fail("no such timed case", argv[i]);
}
}
ph_context_t *ctx = NULL;
if (load_image(&ctx, argv[0])) {
return 1;
}
json_t j = json_begin(stdout);
json_int(&j, "width", ctx->image.width);
json_int(&j, "height", ctx->image.height);
json_string(&j, "build_info", ph_get_build_info());
#ifdef __VERSION__
json_string(&j, "compiler", __VERSION__);
#endif
json_t cases = json_object(&j, "cases");
static double ms[TIME_MAX_RUNS];
for (size_t c = 0; c < COUNT(time_cases); c++) {
const time_case_t *tc = &time_cases[c];
if (!time_chosen(tc->name, argc - 1, argv + 1)) {
continue;
}
int runs = 0;
double total = 0.0;
for (int k = -1; k < TIME_MAX_RUNS; k++) { /* k = -1: the warm-up */
ph_drop_gray_cache(ctx);
double t0 = seconds();
if (tc->run(ctx, argv[0], tc->arg)) {
ph_free(ctx);
return fail("timed case failed", tc->name);
}
double t = seconds() - t0;
if (k >= 0) {
ms[runs++] = t * 1e3;
total += t;
if (runs >= TIME_MIN_RUNS && total >= TIME_MIN_SECONDS) {
break;
}
}
}
qsort(ms, (size_t)runs, sizeof(ms[0]), cmp_double);
json_t o = json_object(&cases, tc->name);
json_double(&o, "min_ms", ms[0]);
json_double(&o, "median_ms", ms[runs / 2]);
json_int(&o, "runs", runs);
json_close_object(&o);
}
json_close_object(&cases);
json_end(&j);
ph_free(ctx);
return 0;
}
The times are measured again only when the library, the tool or the machine changes:
def measure_timing(tool):
"""{"machine": {...}, "small": {...}, "large": {...}, "scan": {...}}, each image's
entry as `site_stages time` prints it (width, height, build_info, compiler, and per
case the minimum and median in milliseconds) and the comparisons of a linear search as
`site_stages scan` prints them (per case, nanoseconds per comparison), from the cache
when nothing that decides them has changed, the machine included."""
host = machine()
paths = [os.path.join(ROOT, p) for p in IMAGES.values()]
key = fingerprint(KEY_FILES, [tool] + paths, json.dumps(host, sort_keys=True))
def measure():
data = {"machine": host}
for name, path in zip(IMAGES, paths):
out = subprocess.run([tool, "time", path], check=True, capture_output=True,
text=True).stdout
data[name] = json.loads(out)
data["scan"] = json.loads(subprocess.run([tool, "scan"], check=True,
capture_output=True, text=True).stdout)
return data
return cached("timing", key, measure, "timing")
Parameters¶
ph_context_set_mhash_params() sets three
values.
| Parameter | Range | Default | Meaning |
|---|---|---|---|
alpha |
above 1 | 2 | the base of the kernel's scale |
level |
0 or above | 1 | its exponent: the kernel reaches \(\lfloor 4\alpha^{\ell} \rfloor\) pixels from its center, at most 32 (a 65×65 kernel) |
size |
62 to 4096 | 512 | the side of step 3; the grid stays 31×31, so a block is \(\lfloor \text{size}/31 \rfloor\) pixels |
The digest is 576 bits whatever the values, and digests are comparable only when computed with the same ones. What the hash sees is decided by one ratio: the kernel's scale against the block's size, both in pixels of the normalized image.
The kernel's scale¶
Raising alpha or level widens the kernel, which then answers to coarser detail. Here
are levels 0, 1 and 2 at alpha 2, their kernels drawn to the same scale, and the block
sums each gives:
The block sums keep their overall pattern, the cactus against the plain wall, and change in detail; a few dozen bits move each way.
The size¶
size changes the block, not the kernel: at 256 a block is 8 pixels, at 1024 it is 33,
against a kernel that stays 17 pixels wide. The coarser the kernel against the block, the
more the block sums follow the picture's broad layout; the finer, the more they follow
texture. At 1024 the digest of the example is nearly as far from the default's as an
unrelated image's.
The work that depends on size grows with its square; the blur, which does not, is the
larger part on a large image. Here is the time against size on the two timed images:
On the small image the time grows slowly up to the default and roughly with the square
of size above it; on the large one the blur's cost hides it up to a size several times
the default.
The numbers behind the chart
| size | 400×400 | 5472×3648 |
|---|---|---|
| 62 | 0.85 ms | 92 ms |
| 128 | 0.91 ms | 94 ms |
| 256 | 1.1 ms | 93 ms |
| 512 | 1.8 ms | 92 ms |
| 1024 | 4.7 ms | 97 ms |
| 2048 | 17 ms | 109 ms |
| 4096 | 66 ms | 159 ms |
The defaults are pHash's, which fixes the size at 512 and does not expose it. Across 24
settings of level and size the property tests' separability moves within a narrow
range with no trend in either, and no setting beats the source's
(provenance has
the sweep).
Settings that affect it¶
Five context settings change the pixels mHash reads.
| Setting | Effect on mHash |
|---|---|
ph_context_set_gray_weights() |
the grayscale formula of step 1 |
ph_context_set_load_grayscale() |
a JPEG decoder converts to grayscale itself (measured), which moves gray levels by one here and there |
ph_context_set_auto_orient() |
on by default: the hash describes the image as displayed, after its EXIF rotation |
ph_context_set_alpha_mode() |
how transparent pixels are composited before grayscale (which background) |
ph_context_set_decode_scale() |
a JPEG read by libjpeg-turbo is decoded at ½, ¼ or ⅛ of its size (measured); the blur's σ is fixed in pixels, so a smaller decode changes what mHash reads, and the digest moves at every scale, toward its threshold for copies |
A gray level here and there reaches the response through the kernel, and the block sums of fine detail are sensitive to it:
On the example photograph the two digests are 14 of 576 bits apart.
In code¶
mHash returns a digest, like BMH, Radial and the two color hashes. This example computes
all five for two images and compares each pair with the function its kind calls for;
mHash's goes through
ph_hamming_distance_digest():
/* Hashes two images with every digest algorithm, compares each pair with the function
* that digest's kind calls for, shows what a mismatched comparison returns, and stores a
* digest as text and reads it back. Build/run:
* cc digest_and_metrics.c -o digest_and_metrics $(pkg-config --cflags --libs libphash)
* ./digest_and_metrics a.jpg b.jpg
*/
#include <libphash.h>
#include <stdio.h>
#include <string.h>
/* Loads `path` into `ctx` in color, so the color hashes can run on it too. */
static int load(ph_context_t *ctx, const char *path) {
ph_error_t err = ph_load_from_file(ctx, path);
if (err != PH_SUCCESS) {
fprintf(stderr, "failed to load %s: %s (%s)\n", path, ph_get_error_string(err),
ph_get_last_error_message(ctx));
return 0;
}
return 1;
}
/* Compares two digests of one algorithm with the function its kind calls for. */
static void compare(const char *name, const ph_digest_t *a, const ph_digest_t *b) {
double score = 0.0;
switch ((ph_digest_kind_t)a->kind) {
case PH_DIGEST_KIND_BITS:
printf("%-13s %3d of %d bits differ, similarity %.3f\n", name,
ph_hamming_distance_digest(a, b), a->size * 8, ph_similarity_digest(a, b));
break;
case PH_DIGEST_KIND_COEFFICIENTS:
if (ph_radial_similarity(a, b, &score) == PH_SUCCESS) {
printf("%-13s peak correlation %.3f (same image at %.1f and above)\n", name, score,
PH_RADIAL_PCC_THRESHOLD);
} else {
printf("%-13s no angular structure to compare\n", name);
}
break;
case PH_DIGEST_KIND_HISTOGRAM:
if (ph_histogram_intersection(a, b, &score) == PH_SUCCESS) {
printf("%-13s histogram intersection %.3f\n", name, score);
}
break;
case PH_DIGEST_KIND_VECTOR16:
case PH_DIGEST_KIND_VECTOR:
printf("%-13s L2 distance %.2f\n", name, ph_l2_distance(a, b));
break;
default:
printf("%-13s kind %d: no comparison\n", name, a->kind);
break;
}
}
int main(int argc, char **argv) {
if (argc != 3) {
fprintf(stderr, "usage: %s <image-a> <image-b>\n", argv[0]);
return 1;
}
static const ph_algorithm_t algos[] = {PH_ALGO_MHASH, PH_ALGO_BMH, PH_ALGO_RADIAL,
PH_ALGO_COLOR_HASH, PH_ALGO_COLOR_MOMENTS};
enum {
N_ALGOS = sizeof(algos) / sizeof(algos[0]),
};
ph_digest_t da[N_ALGOS], db[N_ALGOS];
ph_context_t *ctx = NULL;
if (ph_create(&ctx) != PH_SUCCESS) {
fprintf(stderr, "ph_create failed\n");
return 1;
}
/* One context, reused: each load replaces the previous image and its cached work. */
for (int img = 0; img < 2; img++) {
if (!load(ctx, argv[1 + img])) {
ph_free(ctx);
return 1;
}
for (int i = 0; i < N_ALGOS; i++) {
ph_error_t err = ph_compute_digest(ctx, algos[i], img == 0 ? &da[i] : &db[i]);
if (err != PH_SUCCESS) {
fprintf(stderr, "%s failed: %s\n", ph_algorithm_name(algos[i]),
ph_get_error_string(err));
ph_free(ctx);
return 1;
}
}
}
for (int i = 0; i < N_ALGOS; i++) {
compare(ph_algorithm_name(algos[i]), &da[i], &db[i]);
}
/* The wrong metric is refused, not answered with a meaningless number: a radial
* digest is quantized coefficients, not a bit vector. */
printf("hamming distance of two radial digests: %d (refused)\n",
ph_hamming_distance_digest(&da[2], &db[2]));
/* Store a digest as text and read it back; the text carries the kind. */
char hex[PH_DIGEST_HEX_BUFFER_SIZE];
ph_digest_t back;
if (ph_digest_to_hex(&da[1], hex, sizeof(hex)) != PH_SUCCESS ||
ph_digest_from_hex(hex, &back) != PH_SUCCESS) {
fprintf(stderr, "hex round trip failed\n");
ph_free(ctx);
return 1;
}
printf("BMH as text: %s\n", hex);
printf("read back: %s\n", (back.size == da[1].size && back.kind == da[1].kind &&
memcmp(back.data, da[1].data, back.size) == 0)
? "identical"
: "DIFFERENT");
ph_free(ctx);
return 0;
}
Where it comes from¶
Two sources cover two halves of mHash. The operator is the Laplacian of Gaussian of
D. Marr and E. Hildreth, "Theory of edge detection" (1980), sampled as the Mexican hat.
The hash built on it has no paper: it is ph_mh_imagehash() of the pHash library, by
Evan Klinger and David Starkweather, and Zauner's thesis (2010), which describes it, says
the construction "has not been proposed previously". pHash's code is therefore the primary
source for every step and constant on this page: the blur at 1, the size 512, the 256
levels, the scale \(\alpha = 2\), \(\ell = 1\), the 16-pixel blocks, the 3×3 windows at
stride 4.
Despite the name, neither pHash nor this library detects edges as Marr and Hildreth do, by the zero-crossings of the response: the response is summed over blocks and compared within windows. The operator is theirs; the edge detector is not implemented.
The digest is not bit for bit pHash's. pHash blurs with Deriche's recursive filter and resamples with a quintic filter, where this library uses a direct Gaussian and the Mitchell filter; it also computes the block sums the long way. The construction and every parameter are reproduced; the arithmetic is not (provenance § 5).
-
Measured when this site was built, on Intel(R) Xeon(R) Platinum 8370C CPU @ 2.80GHz, 4 cores, Linux x86_64; compiled by 13.3.0, optimized, with libjpeg-turbo decoding the JPEGs: the minimum of at least five runs after a warm-up. What each time covers, and the code, are under “How this was measured” beside it. ↩