Central difference approximations to estimate a Jacobian matrix, HELP!!!

I have two function F and G defined in script files FM.m and GM.m and I've been asked to use central difference approximations to estimate the Jacobian matrix J =[Fx(1,2) Fy(1,2);Gx(1,2) Gy(1,2)] in which Fx denotes partial differentiation variable.

The central difference approximation is f'(x)=(f(x+h)-f(x-h))/2h.

Any help would be really appreciated, thanks!!