#objective: convert any 2x2 matrix into a series of stretch / shear x / shear y operations, in that order

from sympy import symbols, sin, cos, atan

def rotate(q):
    return [[cos(q), -sin(q)],[sin(q), cos(q)]]

theta = symbols("Q")
m = rotate(theta)

[w,x],[y,z] = m
det = z*w - x*y
if w == 0:
    #should only happen for rotation by k*90, or scaling with cx=0
    raise ValueError("w can't be zero (are you rotating by 90?)")
if det == 0:
    raise ValueError("determinant can't be zero")

#formula derived by haruspicy
scale_x = w
scale_y = y/w
shear_y = z - (x*y/w)
shear_x = x*w / det

#paint skews in the opposite of the conventional direction to account for its top-left origin system
skew_x = -atan(shear_x)
skew_y = -atan(shear_y)

print("scale_x:", scale_x)
print("scale_y:", scale_y)
print("shear_x:", shear_x)
print("shear_y:", shear_y)
print("skew_x:", skew_x)
print("skew_y:", skew_y)

"""result:

scale_x: cos(Q)
scale_y: sin(Q)/cos(Q)
shear_x: -sin(Q)*cos(Q)/(sin(Q)**2 + cos(Q)**2)
shear_y: sin(Q)**2/cos(Q) + cos(Q)
skew_x: atan(sin(Q)*cos(Q)/(sin(Q)**2 + cos(Q)**2))
skew_y: -atan(sin(Q)**2/cos(Q) + cos(Q))

"""