Matrici di trasformazione

26 gennaio 2025

Nell’articolo precedente abbiamo visto come comporre più trasformazioni applicandole in sequenza:

const star = array<vec2f, 30>(
13 collapsed lines
  // Outer star points
  vec2f(0.000000, -0.197140), vec2f(0.309011, -0.425297), vec2f(0.187467, -0.060929),
  vec2f(0.187467, -0.060929), vec2f(0.500000, 0.162443), vec2f(0.115866, 0.159500),
  vec2f(0.115866, 0.159500), vec2f(0.000000, 0.525707), vec2f(-0.115866, 0.159500),
  vec2f(-0.115866, 0.159500), vec2f(-0.500000, 0.162443), vec2f(-0.187467, -0.060929),
  vec2f(-0.187467, -0.060929), vec2f(-0.309011, -0.425297), vec2f(0.000000, -0.197140),

  // Inner star points
  vec2f(0.000000, -0.197140), vec2f(0.000000, 0.000000), vec2f(0.187467, -0.060929),
  vec2f(0.187467, -0.060929), vec2f(0.000000, 0.000000), vec2f(0.115866, 0.159500),
  vec2f(0.115866, 0.159500), vec2f(0.000000, 0.000000), vec2f(-0.115866, 0.159500),
  vec2f(-0.115866, 0.159500), vec2f(0.000000, 0.000000), vec2f(-0.187467, -0.060929),
  vec2f(-0.187467, -0.060929), vec2f(0.000000, 0.000000), vec2f(0.000000, -0.197140)
);

@vertex
fn vs(@builtin(vertex_index) index: u32) -> @builtin(position) vec4f {
  // Step 1: Scala
  let sX = 0.5;
  let sY = 0.5;
  var transformed = vec2f(
    star[index].x * sX,
    star[index].y * sY
  );

  // Step 2: Rotazione
  let theta = radians(45.0);
  let cosT = cos(theta);
  let sinT = sin(theta);
  transformed = vec2f(
    transformed.x * cosT - transformed.y * sinT,
    transformed.x * sinT + transformed.y * cosT
  );

  // Step 3: Traslazione
  let tX = 0.2;
  let tY = -0.1;
  transformed += vec2f(tX, tY);

  return vec4f(transformed, 0.0, 1.0);
}

L’implementazione tuttavia è poco flessibile perchè limitata alla composizione di scala, rotazione e traslazione. Se volessimo applicare delle trasformazioni diverse, o cambiarne l’ordine, saremmo costretti a modificare ogni volta il codice sorgente.

In questo articolo definiremo le trasformazioni geometriche usando le matrici, uno strumento dell’algebra lineare. Potendo rappresentare uniformemente trasformazioni e composizioni, le matrici rendono il codice generico e riutilizzabile per qualunque sequenza di trasformazioni.

Il codice finale è disponibile nel ramo 04-matrix del repository GitHub.

#Matrici

Una matrice rettangolare $M = (m_{ij})$ di dimensioni $m \times n$ è una “tabella” di numeri reali:

M = \begin{bmatrix} m_{11} & m_{12} & \cdots & m_{1n} \\ m_{21} & m_{22} & \cdots & m_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ m_{m1} & m_{m2} & \cdots & m_{mn} \end{bmatrix}

dove $m$ è il numero di righe, $n$ il numero di colonne e $m_{ij}$ sono gli elementi della matrice.

Una matrice di dimensioni $1 \times n$ è detta matrice riga:

\begin{bmatrix} m_{11} & m_{12} & \cdots & m_{1n} \end{bmatrix}

Mentre una matrice di dimensioni $m \times 1$ è detta matrice colonna:

\begin{bmatrix} m_{11} \\ m_{21} \\ \vdots \\ m_{m1} \end{bmatrix}

Infine, una matrice di dimensioni $n \times n$ è detta matrice quadrata di ordine $n$ :

\begin{bmatrix} m_{11} & m_{12} & m_{13} & \cdots & m_{1n} \\ m_{21} & m_{22} & m_{23} & \cdots & m_{2n} \\ m_{31} & m_{32} & m_{33} & \cdots & m_{3n} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ m_{n1} & m_{n2} & m_{n3} & \cdots & m_{nn} \end{bmatrix}

Di seguito sono mostrate in ordine una matrice rettangolare $2 \times 3$ , una matrice quadrata di ordine $2$ e una matrice riga:

\begin{bmatrix} \pi & -e & 1 \\ 0 & 5 & -\frac{1}{2} \end{bmatrix} , \begin{bmatrix} 1 & 2 \\ 3 & 4 \end{bmatrix} , \begin{bmatrix} -1 & 0 & 5 \end{bmatrix}

#Prodotto matriciale

Data una matrice $A = (a_{ij})$ di dimensioni $m \times n$ e una matrice $B = (b_{ij})$ di dimensioni $n \times p$ , chiamiamo prodotto tra $A$ e $B$ una matrice $C = (c_{ij})$ di dimensioni $m \times p$ così definita:

c_{ij} = \sum_{k=1}^{n} a_{ik} b_{kj}

In altre parole, l’elemento $c_{ij}$ si ottiene combinando l’ $i$ -esima riga della matrice $A$ con la $j$ -esima colonna della matrice $B$ . Ad esempio:

\begin{bmatrix} a_{11} & a_{12} \\ a_{21} & a_{22} \end{bmatrix} \begin{bmatrix} b_{11} & b_{12} \\ b_{21} & b_{22} \end{bmatrix} = \begin{bmatrix} a_{11}b_{11} + a_{12}b_{21} & a_{11}b_{12} + a_{12}b_{22} \\ a_{21}b_{11} + a_{22}b_{21} & a_{21}b_{12} + a_{22}b_{22} \end{bmatrix}

Prodotto tra due matrici 2x2. Tocca o clicca sugli elementi per evidenziare le righe/colonne coinvolte nel calcolo.

Dove i singoli elementi sono:

\begin{align*} c_{11} &= a_{11} b_{11} + a_{12} b_{21} \\ c_{12} &= a_{11} b_{12} + a_{12} b_{22} \\ c_{21} &= a_{21} b_{11} + a_{22} b_{21} \\ c_{22} &= a_{21} b_{12} + a_{22} b_{22} \end{align*}

Il prodotto è definito solo tra matrici “compatibili”, cioè quando il numero di colonne della prima matrice coincide con il numero di righe della seconda. Il risultato avrà il numero di righe della prima matrice e di colonne della seconda. Un modo più semplice per ricordare questa regola è l’esempio:

[m \times n] \cdot [n \times p] = [m \times p]

#Trasposizione

Chiamiamo matrice trasposta o semplicemente trasposta di $M$ la matrice $M^T = (m_{ji})$ di dimensioni $n \times m$ ottenuta scambiando le righe con le colonne. Ad esempio:

\begin{pmatrix} 1 & 2 \\ 3 & 4 \\ 5 & 6 \end{pmatrix}^T = \begin{pmatrix} 1 & 3 & 5 \\ 2 & 4 & 6 \end{pmatrix}

La trasposizione del prodotto tra due matrici equivale a trasporre le singole matrici e a scambiarne l’ordine:

(A \cdot B)^T = B^T \cdot A^T

Una matrice $M$ di dimensioni $n \times n$ si dice simmetrica se $M = M^T$ , cioè se coincide con la sua trasposta:

\begin{bmatrix} 3 & 4 \\ 4 & 5 \end{bmatrix}^T = \begin{bmatrix} 3 & 4 \\ 4 & 5 \end{bmatrix}

Gli elementi sono simmetrici rispetto alla diagonale principale ( $\searrow$ ).

#Matrice identità

La matrice di identità è una matrice quadrata in cui tutti gli elementi della diagonale principale valgono $1$ , mentre gli altri $0$ :

I_n = \begin{bmatrix} 1 & 0 & 0 & \cdots & 0 \\ 0 & 1 & 0 & \cdots & 0 \\ 0 & 0 & 1 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & 1 \end{bmatrix}

La matrice d’identità è l’elemento neutro rispetto all’operazione di prodotto, cioè:

A \cdot I_n = I_n \cdot A = A

dove $A$ è una matrice quadrata di ordine $n$ .

#Matrice inversa

Una matrice quadrata di ordine $n$ si dice invertibile se esiste una matrice quadrata $M^{-1}$ di ordine $n$ tale che:

M^{-1} \cdot M = M \cdot M^{-1} = I_n

In tal caso $M^{-1}$ è detta inversa di $M$ . Una matrice non invertibile è anche detta singolare.

#Rappresentazione di un punto

Nelle sezioni successive useremo le matrici per definire le funzioni di trasformazione. Affinché queste trasformazioni, espresse in forma matriciale, possano agire sui punti, è necessario rappresentare i punti stessi in un formato matematico compatibile con le operazioni tra matrici.

Un punto è rappresentato dalle sue coordinate, che possiamo indicare usando un vettore, una matrice riga o una matrice colonna:

(x, y) , \begin{bmatrix} x & y \end{bmatrix} , \begin{bmatrix} x \\ y \end{bmatrix}

Tutte e tre le scritture definiscono lo stesso identico concetto, la differenza è solo a livello di notazione. Pertanto, considereremo equivalenti le tre forme.

#Trasformazioni lineari

Una trasformazione lineare è una funzione nella forma:

f(x, y) = (a x + b y, c x + d y)

dove $a, b, c, d$ sono coefficienti reali. Il nome deriva dal fatto che le coordinate trasformate sono una combinazione lineare delle coordinate del punto originale. Le trasformazioni lineari possono essere espresse come prodotto tra una matrice quadrata di ordine $2$ e una matrice colonna:

f(x, y) = \begin{bmatrix} a & b \\ c & d \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix}

In alternativa, si può definire come prodotto tra una matrice riga e una matrice quadrata di ordine $2$ :

f(x, y) = \begin{bmatrix} x & y \end{bmatrix} \begin{bmatrix} a & c \\ b & d \end{bmatrix}

Nelle sezioni successive descriveremo i punti tramite matrici colonna, quindi useremo la prima forma.

Le trasformazioni lineari preservano il parallelismo tra rette, ma non gli angoli e le aree: ciò significa che rette parallele prima della trasformazione rimarranno parallele anche dopo, ma angoli formati da coppie di rette potrebbero cambiare la loro ampiezza e figure geometriche potrebbero subire una variazione della loro area.

Un quadrato trasformato mediante scala, taglio, rotazione e traslazione. Il parallelismo tra i lati è mantenuto, ma non gli angoli e l'area. La figura tratteggiata è quella originale.

In altre parole, se prendiamo un parallelogramma, dopo una trasformazione lineare esso rimarrà un parallelogramma, ma potrebbe non essere più un rettangolo o un quadrato. Le direzioni delle rette vengono modificate, e le distanze vengono ridimensionate in modo non uniforme, causando la distorsione degli angoli e la variazione delle aree.

#Scala

La trasformazione di scala

f_S(x, y) = (x \cdot s_x, y \cdot s_y)

può essere scritta come:

f_S(x, y) = \begin{bmatrix} s_x & 0 \\ 0 & s_y \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix}

Eseguendo il prodotto matriciale risulta:

\begin{bmatrix} s_x & 0 \\ 0 & s_y \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} = \begin{bmatrix} s_x \cdot x + 0 \cdot y \\ 0 \cdot x + s_y \cdot y \end{bmatrix} = \begin{bmatrix} s_x \cdot x \\ s_y \cdot y \end{bmatrix}

#Rotazione

La rotazione

f_R(x, y) = (x \cdot \cos (\theta) - y \cdot \sin (\theta), x \cdot \sin (\theta) + y \cdot \cos (\theta))

può essere scritta come:

f_R(x, y) = \begin{bmatrix} \cos(\theta) & -\sin(\theta) \\ \sin(\theta) & \cos(\theta) \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix}

Eseguendo il prodotto matriciale risulta:

\begin{bmatrix} \cos(\theta) & -\sin(\theta) \\ \sin(\theta) & \cos(\theta) \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} = \begin{bmatrix} \cos(\theta) \cdot x - \sin(\theta) \cdot y \\ \sin(\theta) \cdot x + \cos(\theta) \cdot y \end{bmatrix}

#Taglio

La trasformazione di taglio

f_{SH}(x, y) = (x + sh_x \cdot y, y + sh_y \cdot x)

può essere scritta come:

f_{SH}(x, y) = \begin{bmatrix} 1 & sh_x \\ sh_y & 1 \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} \\

Eseguendo il prodotto matriciale risulta:

\begin{bmatrix} 1 & sh_x \\ sh_y & 1 \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} = \begin{bmatrix} 1 \cdot x + sh_x \cdot y \\ sh_y \cdot x + 1 \cdot y \end{bmatrix} = \begin{bmatrix} x + sh_x \cdot y \\ y + sh_y \cdot x \end{bmatrix}

#Matrici in WGSL

Il linguaggio WGSL supporta le matrici di dimensioni comprese tra $2 \times 2$ e $4 \times 4$ con i tipi primitivi matCxRf, dove $C$ è il numero di colonne ed $R$ il numero di righe. La notazione è invertita rispetto a quella matematica, dove si indica prima il numero di righe e poi il numero di colonne. Ad esempio, il tipo mat2x4f rappresenta una matrice $4 \times 2$ (quattro righe e due colonne):

\begin{bmatrix} m_{11} & m_{12} \\ m_{21} & m_{22} \\ m_{31} & m_{32} \\ m_{41} & m_{42} \end{bmatrix}

La codifica binaria dei tipi matCxRf è per colonna (in inglese column-major): ciò significa che matCxRf è strutturato in memoria come un array di $C$ vettori di lunghezza $R$ . Ad esempio mat2x4f è concettualmente equivalente a array<vec4f, 2>:

[(m_{11}, m_{21}, m_{31}, m_{41}), \, (m_{12}, m_{22}, m_{32}, m_{42})]

Bisogna prestare molta attenzione a non confondere la notazione matematica delle matrici, che è per righe, con la rappresentazione binaria in WGSL, che è per colonne. Per capire meglio questo concetto, esaminiamo una matrice non simmetrica, come la matrice di taglio:

shearMatrix(sh_x, sh_y) = \begin{bmatrix} 1 & sh_x \\ sh_y & 1 \end{bmatrix}

Il costruttore di mat2x2f segue la rappresentazione column-major e quindi accetta gli elementi per colonna e non per riga:

fn shearMatrix(shx: f32, shy: f32) -> mat2x2f {
  return mat2x2f(1, shy, shx, 1);
}

fn shearMatrixTransposed(shx: f32, shy: f32) -> mat2x2f {
  return mat2x2f(
    1, shx,
    shy, 1
  );
}

L’indentazione della funzione shearMatrixTransposed ricalca visualmente la notazione matematica, che è per righe, ma il risultato in memoria sarà la matrice trasposta:

shearMatrixTransposed(sh_x, sh_y) = \begin{bmatrix} 1 & sh_y \\ sh_x & 1 \end{bmatrix}

Se avessimo considerato una matrice simmetrica, non ci sarebbe stata alcuna differenza, poiché, per definizione, la trasposta di una matrice simmetrica coincide con la matrice stessa.

#Prodotto in WGSL

Il prodotto in WGSL è definito sia tra matrici che tra matrici e vettori. Possiamo aggiornare il nostro shader di esempio sostituendo le trasformazioni di rotazione e scala con le matrici:

fn scaleMatrix(sX: f32, sY: f32) -> mat2x2f {
  return mat2x2f(
    sX, 0,
    0, sY
  );
}

fn rotationMatrix(theta: f32) -> mat2x2f {
  let cosT = cos(theta);
  let sinT = sin(theta);
  return mat2x2f(
    cosT, sinT,
    -sinT, cosT
  );
}

@vertex
fn vs(@builtin(vertex_index) index: u32) -> @builtin(position) vec4f {
  // Step 1: Scala
  var transformed = scaleMatrix(0.5, 0.5) * star[index];

  // Step 2: Rotazione
  transformed = rotationMatrix(radians(45.0)) * transformed;

  // Step 3: Traslazione
  let tX = 0.2;
  let tY = -0.1;
  transformed += vec2f(tX, tY);

  return vec4f(transformed, 0.0, 1.0);
}

Possiamo notare un primo vantaggio dell’uso delle matrici: il codice per applicare due trasformazioni diverse, scala e rotazione, è sempre un prodotto matriciale.

Le funzioni scaleMatrix e rotationMatrix costruiscono le matrici da applicare nella forma $f(x, y) = M \cdot p$ , dove $p$ è un vettore colonna. Il codice dello shader segue questa convenzione e applica i prodotti nella forma M * p anziché p * M.

In generale, in WGSL i vettori non sono considerati nè righe nè colonne, tuttavia entrambe le espressioni M * p e p * M sono valide. Affinché il prodotto sia sempre ben definito, i vettori vengono interpretati di volta in volta come righe o colonne in base all’ordine in cui compaiono nel prodotto con una matrice. Le due espressioni infatti vengono calcolate nel seguente modo:

\begin{align*} M \cdot p = \begin{bmatrix}a & b \\ c & d \end{bmatrix} \begin{bmatrix}x \\ y\end{bmatrix} = \begin{bmatrix} ax + by \\ cx + dy \end{bmatrix} \\ p \cdot M = \begin{bmatrix}x & y\end{bmatrix} \begin{bmatrix}a & b \\ c & d \end{bmatrix} = \begin{bmatrix} ax + cy & bx + dy \end{bmatrix} \end{align*}

I due risultati non differiscono solo per il fatto di essere un vettore riga piuttosto che colonna (irrilevante in WGSL), ma gli elementi in sè sono differenti! Per questo motivo è fondamentale stabilire una convenzione e applicarla sempre allo stesso modo.

Come per l’esempio di codice mostrato, nelle sezioni che seguono adotteremo sempre le matrici da usare nella forma M * p, facendo attenzione al layout column-major nei costruttori matCxRf.

#Composizione di trasformazioni

La versione aggiornata dello shader implementa la formula:

R \cdot (S \cdot p) + t

dove $R$ è la matrice di rotazione, $S$ è la matrice di scala, $p$ è il vertice da trasformare e $t$ è il vettore di traslazione. La sequenza logica dei passaggi è:

Il punto $p$ viene scalato tramite prodotto matriciale con $S$ , ottenendo un primo punto intermedio.
Il punto appena calcolato viene ruotato tramite prodotto matriciale con $R$ , ottenendo così il secondo punto intermedio.
Infine, il secondo punto intermedio viene traslato tramite somma vettoriale con $t$ , ottenendo così il risultato finale.

In particolare, le trasformazioni $R$ ed $S$ vengono eseguite in due passaggi distinti; vediamo ora com’è possibile unificarli. In generale, il prodotto tra matrici gode della proprietà associativa, ovvero:

A \cdot B \cdot C = A \cdot (B \cdot C) = (A \cdot B) \cdot C

Sfruttando tale proprietà nella formula dello shader ( $p$ è una matrice colonna), possiamo scrivere:

(R \cdot S) \cdot p + t

Posto $W = R \cdot S$ possiamo semplificare in:

W \cdot p + t

Il risultato è significativo: la singola matrice $W$ è equivalente alle trasformazioni $S$ ed $R$ in sequenza, cioè ne rappresenta la composizione. A questo punto possiamo semplificare ulteriormente il codice dello shader:

fn scaleMatrix(sX: f32, sY: f32) -> mat2x2f {
4 collapsed lines
  return mat2x2f(
    sX, 0,
    0, sY
  );
}

fn rotationMatrix(theta: f32) -> mat2x2f {
6 collapsed lines
  let cosT = cos(theta);
  let sinT = sin(theta);
  return mat2x2f(
    cosT, sinT,
    -sinT, cosT
  );
}

@vertex
fn vs(@builtin(vertex_index) index: u32) -> @builtin(position) vec4f {
  var worldTransform = rotationMatrix(radians(45.0)) * scaleMatrix(0.5, 0.5);

  // Step 1: Scala + Rotazione
  var transformed = worldTransform * star[index];

  // Step 2: Traslazione
  let tX = 0.2;
  let tY = -0.1;
  transformed += vec2f(tX, tY);

  return vec4f(transformed, 0.0, 1.0);
}

Il vantaggio di questo approccio può apparire limitato nel nostro esempio, dato che il numero complessivo di operazioni non cambia. Tuttavia, nei casi reali worldTransform è solitamente un parametro dello shader e, di conseguenza, non viene ricalcolato per ogni vertice.

Generalizzando, supponiamo di voler applicare nell’ordine una serie di trasformazioni lineari $M_1, M_2, \ldots, M_n$ ad un punto $p$ :

M_n \cdot \ldots \cdot (M_3 \cdot (M_2 \cdot (M_1 \cdot p)))

Sfruttando ripetutamente la proprietà associativa, possiamo rappresentare la composizione tramite una singola matrice $W = M_n \cdot M_{n-1} \cdot \ldots \cdot M_1$ e utilizzarla direttamente:

W \cdot p

Così facendo, la singola matrice $W$ codifica la sequenza originale, astraendo il vertex shader da tipo, numero e ordine delle trasformazioni da applicare. 😮😱

A differenza del prodotto tra numeri reali, il prodotto tra matrici quadrate generalmente non gode della proprietà commutativa, cioè $A \cdot B \neq B \cdot A$ . L’ordine in cui vengono eseguite le trasformazioni è rilevante, così come avevamo visto intuitivamente nell’articolo precedente. È importante notare che l’ordine in cui le matrici compaiono nella definizione di $W$ è l’inverso dell’ordine con cui vengono “applicate” le trasformazioni. Il prodotto $M_n \cdot M_{n-1} \cdot \ldots \cdot M_1$ è equivalente a trasformare prima con $M_1$ , poi con $M_2$ e così via.

#Conclusioni

In questo articolo abbiamo introdotto le matrici e alcune loro proprietà. Abbiamo visto che scala, rotazione e taglio sono trasformazioni lineari e quindi possono essere espresse nella forma $M \cdot p$ , dove $M$ è una matrice quadrata di ordine $2$ .

Tuttavia, fin’ora nella trattazione abbiamo intenzionalmente lasciato da parte la traslazione. Come vedremo infatti, non è possibile rappresentare una traslazione nel piano con le matrici quadrate di ordine $2$ !

Come risolvere questo problema? La risposta nel prossimo articolo, dedicato a spazi proiettivi e coordinate omogenee.

Indietro Trasformazioni geometriche