Loading and Saving NumPy Data

Module 3: NumPy

Loading and Saving Data: The Ghost Camera

Instead of our regular lesson today, we are embarking on a chilling adventure: ghost hunting!

The Ghost Camera

Last night we went ghost hunting with a very special device: a ghost camera. Ghosts are invisible, so we could not see anything in the dark. But the camera beeps when a ghost is nearby, and every time it beeped, we snapped a picture.

By the end of the night, the camera had met so many ghosts that it got completely spooked… and broke. When we downloaded the pictures, every photo was broken, and each one in a very different way.

Luckily, the camera keeps a log. This is what it recorded:

22:13  BEEP! Ghost nearby ............. photo_1 saved
       ERROR: memory shattered. 9 fragments shuffled inside the photo.
23:02  BEEP! Ghost nearby ............. photo_2 saved
       ERROR: sensor frozen. Flattened, halves swapped, contrast lost.
00:31  BEEP! Ghost nearby ............. photo_3 saved
       ERROR: interference on the sensor rows. Photo cut into 4 strips.
01:45  BEEP! Ghost nearby ............. photo_4 saved
       ERROR: lens cursed. Rows scrambled. Photo cut into 5 strips.
       NOTE: 3 lenses found in the camera bag.
03:33  BEEP! BEEP! BEEP! STRONG GHOST . photo_5 saved (colour)
       NOTE: ghost surrounded by small ORANGE flames.
       ERROR: colour channels separated and shuffled. Very low light.
03:34  CAMERA FAILURE. Too many ghosts.

Who did we catch? We do not know yet: the photos will tell us.

Our mission: load the broken files, use the log to understand what happened to each photo, recover all five pictures, and find the case number to close the investigation.

Setting Up Our Toolkit

The photos are stored online, in the course repository. BASE is the address of the folder that contains them. Two helper functions:

show(image) displays an array as an image: grayscale if it is 2D, colour if it is 3D.
check(photo, key) tells you whether you found the key of a photo, for instance check(1, ...), or the case number, check("case", ...).

import hashlib

import matplotlib.pyplot as plt
import numpy as np

# The photos are stored in the course repository
BASE = "https://raw.githubusercontent.com/daniprec/Python-Course/main/data/ghost_photos"


def show(image):
    """Display an image: a 2D array (grayscale, 0 = black, 1 = white) or a 3D array (colour)."""
    if image.ndim == 2:
        plt.imshow(image, cmap="gray", vmin=0, vmax=1)
    else:
        plt.imshow(np.clip(image, 0, 1))
    plt.axis("off")
    plt.show()


_KEYS = {1: 'acba8dd74b724348', 2: 'f5ca38f748a1d6ea', 3: '9ae2bdd7beedc2e7', 4: 'd4735e3a265e16ee', 5: '865736a1c30a82dc', 'case': '2eb2158e4711ef9f'}


def check(photo, key):
    """Check the key of a photo (1 to 5) or of the case ("case")."""
    digest = hashlib.sha256(str(key).strip().upper().encode()).hexdigest()[:16]
    if digest != _KEYS[photo]:
        print(f"Wrong key for photo {photo}. Keep investigating...")
    elif photo == "case":
        print("Correct case number! The ghost hunt is solved.")
    else:
        print(f"Correct key! Photo {photo} is recovered.")

Where Are the Photos?

Each photo is in its own folder. To get the address (URL) of a file, add the folder and the file name to BASE:

url = f"{BASE}/photo_1/shattered.txt"

Photo	Files
1	`photo_1/shattered.txt`
2	`photo_2/frozen.txt`
3	`photo_3/banded.csv`
4	`photo_4/cursed.csv`, and the lenses `photo_4/lens_1.csv`, `lens_2.csv`, `lens_3.csv`
5	`photo_5/colour.txt`

np.loadtxt() reads a text file straight from its URL:

shattered = np.loadtxt(f"{BASE}/photo_1/shattered.txt")

The first time, NumPy downloads the file and keeps a copy next to your notebook, in a folder called raw.githubusercontent.com. Do not be surprised if it appears!

These are the tools you will need. You already know most of them!

Tool	What it does
`np.loadtxt(url)`	Loads a text file (`.txt`); for `.csv` files, add the separator: `delimiter=","` or `delimiter=";"`
`np.save(path, arr)` / `np.load(path)`	Saves / loads an array in NumPy’s binary format (`.npy`), on your computer
`np.savetxt(path, arr)`	Saves an array as a text file
`np.concatenate(arrays, axis=...)`	Joins arrays along an existing axis (also `np.hstack()`, `np.vstack()`)
`np.split(arr, n, axis=...)`	Cuts an array into `n` equal pieces (also `np.hsplit()`, `np.vsplit()`)
`np.stack(arrays, axis=...)`	Joins arrays along a new axis

Every photo is a $300 \times 300$ array, with values between 0 (black) and 1 (white). When a photo is split into 9 fragments, each fragment is $100 \times 100$, and they are placed like this:

(0, 0) (0, 1) (0, 2)
(1, 0) (1, 1) (1, 2)
(2, 0) (2, 1) (2, 2)

Two things are true for every photo of our camera:

It has a white frame: the 4 pixels on the border are 1.
Just inside the white frame, it has a thin black line: those pixels are 0.

The Five Broken Photos

Every photo was saved in a single file, but the camera cut it into pieces and put them back in the wrong order. So recovering a photo always takes two steps:

Rebuild it: split the file into its pieces, and join them again in the right order. Each photo was cut in a different way!
Undo the damage: each photo also suffered its own transformation, and the camera log tells us which one.

This is what we have for each photo, and how to fix it.

Photo 1: The Shattered Photo

Files: photo_1/shattered.txt. One text file, $300 \times 300$.

What happened: the memory shattered. The photo was cut into its nine $100 \times 100$ fragments, and they were put back together in a random order. The photo has the right size, but it looks like a scrambled puzzle.

We name the nine places of the shattered photo with letters, row by row:

A B C
D E F
G H I

How to fix it:

Cut the shattered photo back into its nine fragments: first into three rows of fragments, then each row into three fragments (np.vsplit() and np.hsplit()). Keep them in a list, in the order A (top left fragment), B (top middle fragment), …, I (bottom right fragment).
Find the right order. The white frame tells you which ones go at the corners and at the edges, and the picture does the rest.
Join them again in the right order: first each row of three fragments side by side, then the three rows one on top of the other.

Key: the letters of the places, in the order of the original photo (row by row), for instance check(1, "ABCDEFGHI").

NumPy Cheatsheet

np.hsplit(arr, n) cuts an array into n equal pieces along columns (horizontally).

np.vsplit(arr, n) cuts an array into n equal pieces along rows (vertically).

np.hstack([arr1, arr2, arr3]) joins arrays side by side (horizontally).

np.vstack([arr1, arr2, arr3]) joins arrays one on top of the other (vertically).

np.rot90(arr) rotates a 2D array 90 degrees counterclockwise. np.rot90(arr, k=2) rotates it twice (180 degrees), and np.rot90(arr, k=-1) rotates it clockwise.

Photo 2: The Frozen Photo

Files: photo_2/frozen.txt. One text file… with one value per line.

What happened: the ghost’s chill froze the sensor. The photo was flattened into a single line of values, and its two halves, left and right, were swapped. The photo also lost almost all its contrast: every value was squeezed into a tiny range. On top of that, it was mirrored along the diagonal (rows and columns swapped).

How to fix it:

Rebuild it: give the photo back its square shape, cut it into its left and right halves, and swap them.
Stretch the values back to the range $[0, 1]$. Remember that every photo has pure black (the line) and pure white (the frame): \[\text{stretched} = \frac{\text{stored} - \min}{\max - \min}\]
Undo the mirror.

Tip: the photo is so faint that show() displays a plain gray square. To see what you are doing in step 1, look at the stretched version: stretching does not care about the position of the pixels, so the two steps can go in any order.

Key: how many times did the ghost reduce the contrast? The original photo had a range of 1 (from 0 to 1); compare it with the range of the stored values.

NumPy Cheatsheet

arr.reshape(rows, cols) gives an array a new shape, with the same values.

np.hsplit(arr, n) cuts an array into n equal pieces along columns (horizontally).

np.hstack([arr1, arr2]) joins arrays side by side (horizontally).

arr.min() and arr.max() return the smallest and the largest value of an array.

arr.T is the transpose of a 2D array: it swaps rows and columns.

Photo 3: The Banded Photo

Files: photo_3/banded.csv. One CSV file, $300 \times 300$. The camera wrote it European-style: the values are separated by semicolons (;) instead of commas.

What happened: the ghost interfered with the sensor, row by row. Each row of pixels got darker, but not all by the same amount. All the pixels of a row were multiplied by the same number between 0 and 1, the gain of that row:

a row with gain 1 stayed as it was;
a row with gain 0.5 became half as bright;
a row with gain 0.2 became very dark.

That is why the photo is covered with dark horizontal bands. Pixel by pixel: \[\text{stored}[i, j] = \text{photo}[i, j] \times \text{gain}[i]\]

Then the camera cut the photo into four vertical strips, and shuffled them. This time there are no numbers: you will need your eyes.

How to fix it:

Rebuild it: cut the photo into its four strips and put them in order. The frame helps: look for the strip with the black line on its left edge, and the one with the black line on its right edge. For the two strips in the middle, try both orders and look at the picture.
Find the gains. In the original photo, the first column is part of the white frame: all its values are 1. In row $i$, that pixel became $1 \times \text{gain}[i] = \text{gain}[i]$. So the first column of the rebuilt photo is the list of gains, one per row.
Undo the interference: divide every row by its own gain. Careful with broadcasting! The gains are a 1D array with shape (300,), and NumPy lines up a 1D array with the columns. Turn the gains into a column, with shape (300, 1), so that each row is divided by its gain.

Key: the number (index) of the first row with the strongest interference, that is, the smallest gain.

NumPy Cheatsheet

np.loadtxt(url, delimiter=";") loads a file whose values are separated by semicolons.

np.hsplit(arr, n) cuts an array into n equal pieces along columns (horizontally).

np.hstack([arr1, arr2, arr3]) joins arrays side by side (horizontally).

arr[:, 0] is the first column of a 2D array, with shape (rows,).

arr.reshape(-1, 1) turns a 1D array into a column, with shape (rows, 1).

np.argmin(arr) returns the position of the smallest value (the first one, if there is a tie).

Photo 4: The Cursed Lens

Files: photo_4/cursed.csv, a $300 \times 300$ photo in CSV format (values separated by commas), and lens_1.csv, lens_2.csv, lens_3.csv, three $300 \times 300$ matrices.

What happened: the ghost cursed the lens of the camera. The photo went through a lens matrix $L$, which scrambled its rows: \[\text{stored} = L \cdot \text{photo}\] We found three lenses in the camera bag, but only one of them was used.

Then the camera cut the stored photo into five horizontal strips of the same size, and shuffled them. With the rows scrambled, the strips are impossible to order by eye… but the camera left a trace: the top-left pixel of every strip stores its number, from 1 to 5, counted from top to bottom.

How to fix it:

Rebuild it: cut the photo into its five strips, sort them by their number (np.argsort() is your friend), and join them again.
Undo the matrix product with the inverse of the lens, $\text{photo} = L^{-1} \cdot \text{stored}$. A matrix has an inverse only if its determinant is not 0: that is how you find the lens that was used.

Key: the number of the lens.

NumPy Cheatsheet

np.loadtxt(url, delimiter=",") loads a file whose values are separated by commas.

np.vsplit(arr, n) cuts an array into n equal pieces along rows (vertically).

arr[0, 0] is the top-left value of a 2D array.

np.argsort(arr) returns the positions that would sort the array: the position of the smallest value first.

np.vstack([arr1, arr2, arr3]) joins arrays one on top of the other (vertically).

np.linalg.det(A) returns the determinant of a square matrix.

np.linalg.inv(A) returns the inverse of a square matrix (it fails if the determinant is 0).

A @ B is the matrix product of A and B.

Photo 5: The Colour Photo

Files: photo_5/colour.txt. One text file, $300 \times 900$: three $300 \times 300$ images side by side, numbered 1 to 3 from left to right.

What happened: the strongest ghost of the night, and our only colour photo. A colour image is a 3D array with shape (height, width, 3): one layer, or channel, for red, green, and blue.

The ghost separated the three channels and placed them side by side, in a shuffled order. On top of that, it was night: every value was squared, so the photo is darker than it should be.

How to fix it:

Rebuild it: cut the file into its three images, and join them along a new axis to get a (300, 300, 3) array.
Undo the square.
Which image is red, which is green, and which is blue? Try the possible orders and use the clue from the camera log: this ghost was surrounded by small orange flames. Careful with blue or green flames: those belong to other spirits (will-o’-the-wisps!), not to our ghost.

Key: the numbers of the red, green, and blue images, in that order, for instance check(5, "123").

NumPy Cheatsheet

np.hsplit(arr, n) cuts an array into n equal pieces along columns (horizontally).

np.stack([red, green, blue], axis=2) joins three 2D arrays along a new last axis: three (300, 300) arrays give one (300, 300, 3) array.

arr[:, :, 0] is the first channel (red) of a colour image.

np.sqrt(arr) returns the square root of every value.

The Case Number

Each recovered photo hides one digit of the case number: the second decimal of the mean brightness of photo[100:250, 100:] (rows 100 to 249, columns 100 to the end). For example, a mean of 0.4371 gives the digit 3. For the colour photo, use only its red channel.

Careful: a photo that is only half fixed gives the wrong digit. You need all five photos fully recovered.

Key: the five digits, from photo 1 to photo 5, for instance check("case", "12345").

Your Mission

Work in groups of three. Each member of the group is in charge of different photos (this table is a suggestion, you can split the task in the way you prefer):

Member	Photos	What you will practise
A	1 and 2	Loading, splitting, sorting, joining, `reshape()`, element-wise arithmetic
B	3 and 4	Delimiters, splitting, sorting, broadcasting, matrix products, determinants, inverses
C	5 and the case number	Splitting, a new axis with `np.stack()`, colour images, aggregations

Use check() to confirm each key. When you finish your part, help the rest of your group: the case number needs all five photos, so it can only be found as a team. The first group to find it closes the investigation!

When a photo is recovered, save the evidence on your computer: store it with np.save() and load it back with np.load(), so that you never have to fix it again. Can you also store it with np.savetxt()? Try it with photo 5.

Why Does This Work?

For the curious: the mathematics behind each fix. Click to open.

Photo 2: why can we undo the freeze exactly?

Besides swapping rows and columns (which only moves pixels around), the ghost applied an affine map to every pixel: $\text{stored} = a \cdot \text{photo} + b$, with $a > 0$. An affine map with $a \neq 0$ is invertible, but we need to know $a$ and $b$.

That is where the frame helps. Every photo contains pure black (0) and pure white (1), so the smallest stored value is $b$ and the largest is $a + b$. The stretch $\frac{\text{stored} - \min}{\max - \min}$ recovers the photo exactly. Without pure black and pure white in the photo, the stretch would make it too contrasted, and the digit would be wrong.

Photo 3: broadcasting is a matrix product in disguise

Multiplying every row $i$ by its gain $g_i$ is the same as multiplying on the left by a diagonal matrix:

\[\text{stored} = G \cdot \text{photo}, \qquad G = \begin{pmatrix} g_1 & & \\ & \ddots & \\ & & g_{300} \end{pmatrix}\]

The inverse of a diagonal matrix is diagonal, with entries $1/g_i$. So stored / gain.reshape(-1, 1) computes exactly $G^{-1} \cdot \text{stored}$, without building a $300 \times 300$ matrix. And stored / gain, the wrong version, computes $\text{stored} \cdot G^{-1}$: it divides the columns. Try it: np.linalg.inv(np.diag(gain)) @ stored.

Photo 4: why is the inverse of the lens its transpose?

The correct lens is a permutation matrix: the identity matrix with its rows shuffled. Each row and each column has a single 1, so the rows are orthonormal vectors. That makes $L$ an orthogonal matrix: $L L^T = I$, and therefore $L^{-1} = L^T$.

Its determinant is $\pm 1$. Each swap of two rows of the identity changes the sign of the determinant, and $|\det L| = 1$ means that the lens does not stretch or squash space: it only moves the rows around.

Photo 4: what does a determinant of 0 mean here?

A matrix with determinant 0 loses information, so no matrix can undo it.

Lens 1 has two equal rows: two rows of the stored photo would be copies of the same row of the original, and one row of the original would be lost forever.
Lens 3 has a row of zeros: one row of the stored photo would be black, whatever the original was.

In both cases, different photos would give the same stored photo, so it is impossible to know which one we started from.

Photo 5: why the square root, and why six orders?

On $[0, 1]$, the function $x \mapsto x^2$ is increasing, so it is invertible, and its inverse is $\sqrt{x}$. Since $x^2 \le x$ there, squaring darkens every pixel, except pure black and pure white.

The three channels can be placed in $3! = 3 \times 2 \times 1 = 6$ orders: a permutation of the last axis. The mean brightness of the photo does not depend on the order, but the mean of the red channel does: that is why the digit of photo 5 uses the red channel only.

What Else Can We Load and Save?

Today we only needed np.loadtxt(), because the photos were tidy text files on the internet. When you work with files on your own computer, NumPy and Python offer much more:

Skip headers and pick columns. Real files often start with a line of column names. skiprows skips lines, and usecols keeps only some columns:
```
data = np.loadtxt("readings.csv", delimiter=",", skiprows=1, usecols=(0, 2))
```
Missing values. np.loadtxt() fails if a value is missing. np.genfromtxt() reads it anyway and writes nan (not a number) in the empty places:
```
data = np.genfromtxt("readings.csv", delimiter=",", skip_header=1)
print(np.isnan(data).sum())   # how many values are missing
```
Readable files for other programs. np.savetxt() can choose the separator, the number format, and a header, so that Excel can open the file:
```
np.savetxt("photo_1.csv", photo_1, delimiter=",", fmt="%.3f", header="my recovered photo")
```
Several arrays in one file. np.savez() stores many arrays together, each with a name; np.savez_compressed() does the same and takes less space:
```
np.savez("evidence.npz", photo_1=photo_1, photo_5=photo_5)
evidence = np.load("evidence.npz")
print(evidence["photo_5"].shape)
```
Images. Matplotlib reads and writes image files as NumPy arrays: plt.imread("ghost.png") and plt.imsave("ghost.png", photo_5). PNG images often have a fourth channel for transparency, so their shape is (height, width, 4).
Many files. On your computer, Python can list the files of a folder for you: os.listdir("ghost_photos/photo_2") or glob.glob("ghost_photos/photo_2/*.txt"). Today we did not need it, because we knew the names of the files.

Summary and Reflection

np.loadtxt() loads text files, even straight from a URL; delimiter sets the separator of CSV files.
np.load() loads .npy files from your computer.
np.save() keeps an array exactly, with its shape and data type; np.savetxt() writes readable text, but only for 1D and 2D arrays.
np.concatenate(), np.hstack(), and np.vstack() join arrays along an existing axis; np.stack() creates a new one. np.split(), np.hsplit(), and np.vsplit() undo the joining.
Array arithmetic, broadcasting, aggregations, and matrix operations can undo many transformations, as long as we know what happened.

Which photo was the hardest to recover, and why? For each photo, how did you know that your result was correct?

Solutions

After you finish your investigation, check your methods against the proposed solutions.

Solutions

Homework and Review

Work through the homework exercises available here.

Homework

--- title: "Loading and Saving NumPy Data" subtitle: "Module 3: NumPy" format: html execute: enabled: false # TEACHER NOTE (not rendered): the ghost photos for this session. # - Source images and generator: scripts/ghost_hunt/ in this repo # (ghost_1.png ... ghost_5.png and generate_ghost_photos.py). # - To regenerate: python scripts/ghost_hunt/generate_ghost_photos.py (needs numpy and Pillow). # It also writes data/ghost_photos.zip, the offline backup linked in the page. # It writes the files of data/ghost_photos/ (read by the students from their GitHub URLs), # plus scripts/ghost_hunt/answers.json with the keys. # - answers.json is in .gitignore: never commit it. The fixed seed always gives the same # photos and keys. --- # Loading and Saving Data: The Ghost Camera Instead of our regular lesson today, we are embarking on a _chilling_ adventure: **ghost hunting!** ## The Ghost Camera Last night we went ghost hunting with a very special device: a **ghost camera**. Ghosts are invisible, so we could not see anything in the dark. But the camera beeps when a ghost is nearby, and every time it beeped, we snapped a picture. By the end of the night, the camera had met so many ghosts that it got completely spooked... and broke. When we downloaded the pictures, every photo was broken, and each one in a **very different way**. Luckily, the camera keeps a log. This is what it recorded: ```text 22:13 BEEP! Ghost nearby ............. photo_1 saved ERROR: memory shattered. 9 fragments shuffled inside the photo. 23:02 BEEP! Ghost nearby ............. photo_2 saved ERROR: sensor frozen. Flattened, halves swapped, contrast lost. 00:31 BEEP! Ghost nearby ............. photo_3 saved ERROR: interference on the sensor rows. Photo cut into 4 strips. 01:45 BEEP! Ghost nearby ............. photo_4 saved ERROR: lens cursed. Rows scrambled. Photo cut into 5 strips. NOTE: 3 lenses found in the camera bag. 03:33 BEEP! BEEP! BEEP! STRONG GHOST . photo_5 saved (colour) NOTE: ghost surrounded by small ORANGE flames. ERROR: colour channels separated and shuffled. Very low light. 03:34 CAMERA FAILURE. Too many ghosts. ``` Who did we catch? We do not know yet: the photos will tell us. Our mission: load the broken files, use the log to understand what happened to each photo, **recover all five pictures**, and find the **case number** to close the investigation. --- ## Setting Up Our Toolkit The photos are stored online, in the course repository. `BASE` is the address of the folder that contains them. Two helper functions: - `show(image)` displays an array as an image: grayscale if it is 2D, colour if it is 3D. - `check(photo, key)` tells you whether you found the key of a photo, for instance `check(1, ...)`, or the case number, `check("case", ...)`. ```{python} import hashlib import matplotlib.pyplot as plt import numpy as np # The photos are stored in the course repository BASE = "https://raw.githubusercontent.com/daniprec/Python-Course/main/data/ghost_photos" def show(image): """Display an image: a 2D array (grayscale, 0 = black, 1 = white) or a 3D array (colour).""" if image.ndim == 2: plt.imshow(image, cmap="gray", vmin=0, vmax=1) else: plt.imshow(np.clip(image, 0, 1)) plt.axis("off") plt.show() _KEYS = {1: 'acba8dd74b724348', 2: 'f5ca38f748a1d6ea', 3: '9ae2bdd7beedc2e7', 4: 'd4735e3a265e16ee', 5: '865736a1c30a82dc', 'case': '2eb2158e4711ef9f'} def check(photo, key): """Check the key of a photo (1 to 5) or of the case ("case").""" digest = hashlib.sha256(str(key).strip().upper().encode()).hexdigest()[:16] if digest != _KEYS[photo]: print(f"Wrong key for photo {photo}. Keep investigating...") elif photo == "case": print("Correct case number! The ghost hunt is solved.") else: print(f"Correct key! Photo {photo} is recovered.") ``` ## Where Are the Photos? Each photo is in its own folder. To get the address (URL) of a file, add the folder and the file name to `BASE`: ```python url = f"{BASE}/photo_1/shattered.txt" ``` | Photo | Files | |-------|-------| | 1 | `photo_1/shattered.txt` | | 2 | `photo_2/frozen.txt` | | 3 | `photo_3/banded.csv` | | 4 | `photo_4/cursed.csv`, and the lenses `photo_4/lens_1.csv`, `lens_2.csv`, `lens_3.csv` | | 5 | `photo_5/colour.txt` | `np.loadtxt()` reads a text file straight from its URL: ```python shattered = np.loadtxt(f"{BASE}/photo_1/shattered.txt") ``` The first time, NumPy downloads the file and keeps a copy next to your notebook, in a folder called `raw.githubusercontent.com`. Do not be surprised if it appears! These are the tools you will need. You already know most of them! | Tool | What it does | |------|--------------| | `np.loadtxt(url)` | Loads a text file (`.txt`); for `.csv` files, add the separator: `delimiter=","` or `delimiter=";"` | | `np.save(path, arr)` / `np.load(path)` | Saves / loads an array in NumPy's binary format (`.npy`), on your computer | | `np.savetxt(path, arr)` | Saves an array as a text file | | `np.concatenate(arrays, axis=...)` | Joins arrays along an existing axis (also `np.hstack()`, `np.vstack()`) | | `np.split(arr, n, axis=...)` | Cuts an array into `n` equal pieces (also `np.hsplit()`, `np.vsplit()`) | | `np.stack(arrays, axis=...)` | Joins arrays along a **new** axis | Every photo is a $300 \times 300$ array, with values between 0 (black) and 1 (white). When a photo is split into 9 fragments, each fragment is $100 \times 100$, and they are placed like this: ```text (0, 0) (0, 1) (0, 2) (1, 0) (1, 1) (1, 2) (2, 0) (2, 1) (2, 2) ``` Two things are true for **every** photo of our camera: - It has a **white frame**: the 4 pixels on the border are `1`. - Just inside the white frame, it has a thin **black line**: those pixels are `0`. --- ## The Five Broken Photos Every photo was saved in **a single file**, but the camera cut it into pieces and put them back in the wrong order. So recovering a photo always takes two steps: 1. **Rebuild it:** split the file into its pieces, and join them again in the right order. Each photo was cut in a different way! 2. **Undo the damage:** each photo also suffered its own transformation, and the camera log tells us which one. This is what we have for each photo, and how to fix it. ### Photo 1: The Shattered Photo **Files:** `photo_1/shattered.txt`. One text file, $300 \times 300$. **What happened:** the memory shattered. The photo was cut into its nine $100 \times 100$ fragments, and they were put back together **in a random order**. The photo has the right size, but it looks like a scrambled puzzle. We name the nine places of the shattered photo with letters, row by row: ```text A B C D E F G H I ``` **How to fix it:** 1. Cut the shattered photo back into its nine fragments: first into three rows of fragments, then each row into three fragments (`np.vsplit()` and `np.hsplit()`). Keep them in a list, in the order A (top left fragment), B (top middle fragment), ..., I (bottom right fragment). 2. Find the right order. The white frame tells you which ones go at the corners and at the edges, and the picture does the rest. 3. Join them again in the right order: first each row of three fragments side by side, then the three rows one on top of the other. **Key:** the letters of the places, in the order of the original photo (row by row), for instance `check(1, "ABCDEFGHI")`. ::: {.callout-tip collapse="true"} ## NumPy Cheatsheet `np.hsplit(arr, n)` cuts an array into `n` equal pieces along **columns** (horizontally). `np.vsplit(arr, n)` cuts an array into `n` equal pieces along **rows** (vertically). `np.hstack([arr1, arr2, arr3])` joins arrays side by side (horizontally). `np.vstack([arr1, arr2, arr3])` joins arrays one on top of the other (vertically). `np.rot90(arr)` rotates a 2D array 90 degrees **counterclockwise**. `np.rot90(arr, k=2)` rotates it twice (180 degrees), and `np.rot90(arr, k=-1)` rotates it clockwise. ::: ### Photo 2: The Frozen Photo **Files:** `photo_2/frozen.txt`. One text file... with **one value per line**. **What happened:** the ghost's chill froze the sensor. The photo was flattened into a single line of values, and its two halves, left and right, were **swapped**. The photo also lost almost all its contrast: every value was squeezed into a tiny range. On top of that, it was **mirrored along the diagonal** (rows and columns swapped). **How to fix it:** 1. Rebuild it: give the photo back its square shape, cut it into its left and right halves, and swap them. 2. Stretch the values back to the range $[0, 1]$. Remember that every photo has pure black (the line) and pure white (the frame): $$\text{stretched} = \frac{\text{stored} - \min}{\max - \min}$$ 3. Undo the mirror. Tip: the photo is so faint that `show()` displays a plain gray square. To see what you are doing in step 1, look at the stretched version: stretching does not care about the position of the pixels, so the two steps can go in any order. **Key:** how many times did the ghost reduce the contrast? The original photo had a range of 1 (from 0 to 1); compare it with the range of the stored values. ::: {.callout-tip collapse="true"} ## NumPy Cheatsheet `arr.reshape(rows, cols)` gives an array a new shape, with the same values. `np.hsplit(arr, n)` cuts an array into `n` equal pieces along **columns** (horizontally). `np.hstack([arr1, arr2])` joins arrays side by side (horizontally). `arr.min()` and `arr.max()` return the smallest and the largest value of an array. `arr.T` is the **transpose** of a 2D array: it swaps rows and columns. ::: ### Photo 3: The Banded Photo **Files:** `photo_3/banded.csv`. One CSV file, $300 \times 300$. The camera wrote it European-style: the values are separated by **semicolons** (`;`) instead of commas. **What happened:** the ghost interfered with the sensor, **row by row**. Each row of pixels got darker, but not all by the same amount. All the pixels of a row were multiplied by the same number between 0 and 1, the **gain** of that row: - a row with gain `1` stayed as it was; - a row with gain `0.5` became half as bright; - a row with gain `0.2` became very dark. That is why the photo is covered with dark horizontal bands. Pixel by pixel: $$\text{stored}[i, j] = \text{photo}[i, j] \times \text{gain}[i]$$ Then the camera cut the photo into **four vertical strips**, and shuffled them. This time there are no numbers: you will need your eyes. **How to fix it:** 1. Rebuild it: cut the photo into its four strips and put them in order. The frame helps: look for the strip with the black line on its left edge, and the one with the black line on its right edge. For the two strips in the middle, try both orders and look at the picture. 2. Find the gains. In the original photo, the first column is part of the white frame: all its values are `1`. In row $i$, that pixel became $1 \times \text{gain}[i] = \text{gain}[i]$. So the first column of the rebuilt photo **is** the list of gains, one per row. 3. Undo the interference: divide every row by its own gain. Careful with broadcasting! The gains are a 1D array with shape `(300,)`, and NumPy lines up a 1D array with the **columns**. Turn the gains into a column, with shape `(300, 1)`, so that each row is divided by its gain. **Key:** the number (index) of the first row with the strongest interference, that is, the smallest gain. ::: {.callout-tip collapse="true"} ## NumPy Cheatsheet `np.loadtxt(url, delimiter=";")` loads a file whose values are separated by semicolons. `np.hsplit(arr, n)` cuts an array into `n` equal pieces along **columns** (horizontally). `np.hstack([arr1, arr2, arr3])` joins arrays side by side (horizontally). `arr[:, 0]` is the first column of a 2D array, with shape `(rows,)`. `arr.reshape(-1, 1)` turns a 1D array into a **column**, with shape `(rows, 1)`. `np.argmin(arr)` returns the **position** of the smallest value (the first one, if there is a tie). ::: ### Photo 4: The Cursed Lens **Files:** `photo_4/cursed.csv`, a $300 \times 300$ photo in CSV format (values separated by commas), and `lens_1.csv`, `lens_2.csv`, `lens_3.csv`, three $300 \times 300$ matrices. **What happened:** the ghost cursed the lens of the camera. The photo went through a lens matrix $L$, which scrambled its rows: $$\text{stored} = L \cdot \text{photo}$$ We found three lenses in the camera bag, but only one of them was used. Then the camera cut the stored photo into **five horizontal strips** of the same size, and shuffled them. With the rows scrambled, the strips are impossible to order by eye... but the camera left a trace: the **top-left pixel** of every strip stores its number, from 1 to 5, counted from top to bottom. **How to fix it:** 1. Rebuild it: cut the photo into its five strips, sort them by their number (`np.argsort()` is your friend), and join them again. 2. Undo the matrix product with the inverse of the lens, $\text{photo} = L^{-1} \cdot \text{stored}$. A matrix has an inverse only if its determinant is not 0: that is how you find the lens that was used. **Key:** the number of the lens. ::: {.callout-tip collapse="true"} ## NumPy Cheatsheet `np.loadtxt(url, delimiter=",")` loads a file whose values are separated by commas. `np.vsplit(arr, n)` cuts an array into `n` equal pieces along **rows** (vertically). `arr[0, 0]` is the top-left value of a 2D array. `np.argsort(arr)` returns the **positions** that would sort the array: the position of the smallest value first. `np.vstack([arr1, arr2, arr3])` joins arrays one on top of the other (vertically). `np.linalg.det(A)` returns the determinant of a square matrix. `np.linalg.inv(A)` returns the inverse of a square matrix (it fails if the determinant is 0). `A @ B` is the matrix product of `A` and `B`. ::: ### Photo 5: The Colour Photo **Files:** `photo_5/colour.txt`. One text file, $300 \times 900$: three $300 \times 300$ images side by side, numbered 1 to 3 from left to right. **What happened:** the strongest ghost of the night, and our only **colour** photo. A colour image is a 3D array with shape `(height, width, 3)`: one layer, or *channel*, for red, green, and blue. The ghost separated the three channels and placed them side by side, in a **shuffled order**. On top of that, it was night: every value was squared, so the photo is darker than it should be. **How to fix it:** 1. Rebuild it: cut the file into its three images, and join them along a **new** axis to get a `(300, 300, 3)` array. 2. Undo the square. 3. Which image is red, which is green, and which is blue? Try the possible orders and use the clue from the camera log: this ghost was surrounded by small **orange flames**. Careful with blue or green flames: those belong to other spirits (will-o'-the-wisps!), not to our ghost. **Key:** the numbers of the red, green, and blue images, in that order, for instance `check(5, "123")`. ::: {.callout-tip collapse="true"} ## NumPy Cheatsheet `np.hsplit(arr, n)` cuts an array into `n` equal pieces along **columns** (horizontally). `np.stack([red, green, blue], axis=2)` joins three 2D arrays along a **new** last axis: three `(300, 300)` arrays give one `(300, 300, 3)` array. `arr[:, :, 0]` is the first channel (red) of a colour image. `np.sqrt(arr)` returns the square root of every value. ::: ### The Case Number Each recovered photo hides one digit of the case number: the **second decimal of the mean brightness of `photo[100:250, 100:]`** (rows 100 to 249, columns 100 to the end). For example, a mean of `0.4371` gives the digit `3`. For the colour photo, use only its red channel. Careful: a photo that is only half fixed gives the wrong digit. You need all five photos fully recovered. **Key:** the five digits, from photo 1 to photo 5, for instance `check("case", "12345")`. --- ## Your Mission Work in **groups of three**. Each member of the group is in charge of different photos (this table is a suggestion, you can split the task in the way you prefer): | Member | Photos | What you will practise | |--------|--------|------------------------| | A | 1 and 2 | Loading, splitting, sorting, joining, `reshape()`, element-wise arithmetic | | B | 3 and 4 | Delimiters, splitting, sorting, broadcasting, matrix products, determinants, inverses | | C | 5 and the case number | Splitting, a new axis with `np.stack()`, colour images, aggregations | Use `check()` to confirm each key. When you finish your part, help the rest of your group: the **case number** needs all five photos, so it can only be found as a team. The first group to find it closes the investigation! When a photo is recovered, **save the evidence** on your computer: store it with `np.save()` and load it back with `np.load()`, so that you never have to fix it again. Can you also store it with `np.savetxt()`? Try it with photo 5. --- ## Why Does This Work? For the curious: the mathematics behind each fix. Click to open. ::: {.callout-tip collapse="true"} ## Photo 2: why can we undo the freeze exactly? Besides swapping rows and columns (which only moves pixels around), the ghost applied an **affine** map to every pixel: $\text{stored} = a \cdot \text{photo} + b$, with $a > 0$. An affine map with $a \neq 0$ is invertible, but we need to know $a$ and $b$. That is where the frame helps. Every photo contains pure black (`0`) and pure white (`1`), so the smallest stored value is $b$ and the largest is $a + b$. The stretch $\frac{\text{stored} - \min}{\max - \min}$ recovers the photo exactly. Without pure black and pure white in the photo, the stretch would make it too contrasted, and the digit would be wrong. ::: ::: {.callout-tip collapse="true"} ## Photo 3: broadcasting is a matrix product in disguise Multiplying every row $i$ by its gain $g_i$ is the same as multiplying on the left by a **diagonal matrix**: $$\text{stored} = G \cdot \text{photo}, \qquad G = \begin{pmatrix} g_1 & & \\ & \ddots & \\ & & g_{300} \end{pmatrix}$$ The inverse of a diagonal matrix is diagonal, with entries $1/g_i$. So `stored / gain.reshape(-1, 1)` computes exactly $G^{-1} \cdot \text{stored}$, without building a $300 \times 300$ matrix. And `stored / gain`, the wrong version, computes $\text{stored} \cdot G^{-1}$: it divides the columns. Try it: `np.linalg.inv(np.diag(gain)) @ stored`. ::: ::: {.callout-tip collapse="true"} ## Photo 4: why is the inverse of the lens its transpose? The correct lens is a **permutation matrix**: the identity matrix with its rows shuffled. Each row and each column has a single `1`, so the rows are orthonormal vectors. That makes $L$ an **orthogonal** matrix: $L L^T = I$, and therefore $L^{-1} = L^T$. Its determinant is $\pm 1$. Each swap of two rows of the identity changes the sign of the determinant, and $|\det L| = 1$ means that the lens does not stretch or squash space: it only moves the rows around. ::: ::: {.callout-tip collapse="true"} ## Photo 4: what does a determinant of 0 mean here? A matrix with determinant 0 loses information, so no matrix can undo it. - Lens 1 has **two equal rows**: two rows of the stored photo would be copies of the same row of the original, and one row of the original would be lost forever. - Lens 3 has a **row of zeros**: one row of the stored photo would be black, whatever the original was. In both cases, different photos would give the same stored photo, so it is impossible to know which one we started from. ::: ::: {.callout-tip collapse="true"} ## Photo 5: why the square root, and why six orders? On $[0, 1]$, the function $x \mapsto x^2$ is increasing, so it is invertible, and its inverse is $\sqrt{x}$. Since $x^2 \le x$ there, squaring darkens every pixel, except pure black and pure white. The three channels can be placed in $3! = 3 \times 2 \times 1 = 6$ orders: a permutation of the last axis. The mean brightness of the photo does not depend on the order, but the mean of the red channel does: that is why the digit of photo 5 uses the red channel only. ::: --- ## What Else Can We Load and Save? Today we only needed `np.loadtxt()`, because the photos were tidy text files on the internet. When you work with files **on your own computer**, NumPy and Python offer much more: - **Skip headers and pick columns.** Real files often start with a line of column names. `skiprows` skips lines, and `usecols` keeps only some columns: ```python data = np.loadtxt("readings.csv", delimiter=",", skiprows=1, usecols=(0, 2)) ``` - **Missing values.** `np.loadtxt()` fails if a value is missing. `np.genfromtxt()` reads it anyway and writes `nan` (not a number) in the empty places: ```python data = np.genfromtxt("readings.csv", delimiter=",", skip_header=1) print(np.isnan(data).sum()) # how many values are missing ``` - **Readable files for other programs.** `np.savetxt()` can choose the separator, the number format, and a header, so that Excel can open the file: ```python np.savetxt("photo_1.csv", photo_1, delimiter=",", fmt="%.3f", header="my recovered photo") ``` - **Several arrays in one file.** `np.savez()` stores many arrays together, each with a name; `np.savez_compressed()` does the same and takes less space: ```python np.savez("evidence.npz", photo_1=photo_1, photo_5=photo_5) evidence = np.load("evidence.npz") print(evidence["photo_5"].shape) ``` - **Images.** Matplotlib reads and writes image files as NumPy arrays: `plt.imread("ghost.png")` and `plt.imsave("ghost.png", photo_5)`. PNG images often have a fourth channel for transparency, so their shape is `(height, width, 4)`. - **Many files.** On your computer, Python can list the files of a folder for you: `os.listdir("ghost_photos/photo_2")` or `glob.glob("ghost_photos/photo_2/*.txt")`. Today we did not need it, because we knew the names of the files. --- ## Summary and Reflection - `np.loadtxt()` loads text files, even straight from a URL; `delimiter` sets the separator of CSV files. - `np.load()` loads `.npy` files from your computer. - `np.save()` keeps an array exactly, with its shape and data type; `np.savetxt()` writes readable text, but only for 1D and 2D arrays. - `np.concatenate()`, `np.hstack()`, and `np.vstack()` join arrays along an existing axis; `np.stack()` creates a new one. `np.split()`, `np.hsplit()`, and `np.vsplit()` undo the joining. - Array arithmetic, broadcasting, aggregations, and matrix operations can undo many transformations, as long as we know what happened. Which photo was the hardest to recover, and why? For each photo, how did you know that your result was correct? --- ## Solutions After you finish your investigation, check your methods against the proposed solutions. [Solutions](load_solutions.qmd#loading-and-saving-data){.btn .btn-primary} --- # Homework and Review Work through the homework exercises available here. [Homework](homework_numpy.qmd#loading-and-saving-data){.btn .btn-primary}