If you were just rotating about the origin, the transformation matrix would just be your standard

But now you have to perform 3 operations:

1) Make (2;1) the "new" origin by shifting it to be the origin

2) Rotate normally about the origin

3) Translate back to the origin of (2;1)

So your resulting transformation matrix will be the product of these three separate transformation matrices. The second is just the standard rotation matrix above (with theta = 60 degrees in your case).

The first transformation, using homogeneous coordinates, translating (2;1) to be the new origin, is

[1 0 -2

0 1 -1

0 0 1]

And the third transformation, translating the origin to (2,1), is

[1 0 2

0 1 1

0 0 1]

So, multiplying these three matrices together will yield the final transformation matrix.