I solved the problem but before I post the solution / explanation, I'll address the previous two replies first (thanks for taking the time and effort to respond).
That's already a bit confusing. You are going to implement transformations in some affine space. So far, so good, but I am not sure what it means to performing them "around the current coordinate system". I'll ignore that for a minute.
I should have been clearer to explain what I mean by this (everything becomes clearer afterwards) :)
Obviously, there's only one coordinate system and we just have point * CTM at the end of the day but visually, a matrix can be thought of as defining another coordinate system. That coordinate system's origin is given by mapping the point (0, 0) using the matrix, the x-axis' unit vector is defined by mapping the points (1, 0) and (0, 0) and the y-axis' unit vector is defined by mapping the points (0, 1) and (0, 0).
So, "scale Y around current coordinate system" means "scale around point given by mapping (0, 0) using the matrix and in the the direction of the y axis unit vector (as defined above)".
That's fine, but it's not what "row major" means.
Sorry, I used the wrong term: I meant points being row vectors.
Performing a -90º rotation and then scaling Y by 2 cannot possibly give you the same result than doing the operations in the reverse order, although I don't understand what you mean by the bit about the CTM.
Obviously, matrix multiplication is not commutative but that's not what we have here due to my definition of scale being performed "around the current coordinate system" (probably not the best wording). If we firstly rotated and then scaled, the scale operation would take into account the rotation matrix and perform the scale in the appropriate direct (e.g., Y scale does not mean actually Y, it means Y in the direction defined by the CTM).
I'm not sure which is which offhand, but I'm pretty sure "inside" is what you want. Then a Y scale matrix will scale the local Y axis, so it will have the same effect whether it's done before or after a rotation.
Correct: inside is what geometrically can be thought of as transforming the coordinate system. But it does not quite work when rotations get involved. The reason is that if you firstly scale, we get "scale * CTM" and then we rotate, the rotation will be prepended, i.e., "rotate * scale * CTM" and the scale no longer has the same meaning when it was applied (e.g., if scale was along Y but rotation was 90º, the scale would instead end up rotating the logical x-axis).
In any case, here are the functions that I've implemented which achieve the desired effect. Each of them takes the CTM to be able to construct the matrix according to the definition.
function scale(ctm, sx, sy) {
scale_x = scale_along_vector(ctm, (1, 0), sx)
scale_y = scale_along_vector(ctm, (0, 1), sy)
return scale_x * scale_y
}
function scale_along_vector(ctm, p, k) {
// k = scale
mapped_origin = (0, 0) * ctm
mapped_point = p * ctm
v = normalise (mapped_point - mapped_origin)
// 3x3 matrix in form of:
// (a b 0)
// (c d 0)
// (0 0 1)
scale_transform = {
.a = 1 + (k - 1) * v.x^2
.b = (k - 1) * v.x * v.y
.c = (k - 1) * v.x * v.y
.d = 1 + (k - 1) * v.y^2
};
translate_to_origin = translation(-mapped_origin.x, -mapped_origin.y)
translate_back = translation(mapped_origin.x, mapped_origin.y)
return translate_to_origin * scale_transform * translate_back
}
function rotate(ctm, radians) {
// determine multiple (1.0 or -1.0) by mapping (1, 0) and
// (0, 1) using the ctm and then determining whether the
// rotation should be in opposite or normal direction
multiplier = ...
rotation_transform = standard rotation matrix(radians * multiplier)
translate_to_origin = translation(-mapped_origin.x, -mapped_origin.y)
translate_back = translation(mapped_origin.x, mapped_origin.y)
return translate_to_origin * rotation_transform * translate_back
}
function translate(ctm, tx, ty) {
// compute unit vectors
mapped_origin = (0, 0) * ctm
xv = ((1, 0) * ctm) - mapped_origin
yv = ((0, 1) * ctm) - mapped_origin
// translate by tx, ty along unit vectors
dx = xv.x * tx + yv.x * ty
dy = xv.y * tx + yv.y * ty
return translation matrix(dx, dy)
}