Wind Stress Curl from wind data

% USAGE: curlz=ra_windstrcurl(lat,lon,u,v)
%
% PREREQUISITE:
% ra_windstr.m written by same author
%
% DESCRIPTION: Function to compute wind stress curl from wind field data Based on
% NRSC (2013)
%
% INPUTS:
% lat = Latitude vector [deg.]
% lon = Longitude Vector [deg.]
% u = Zonal wind component [m/s], must be 2D
% v = Meridional wind component [m/s], must be 2D
%
% OUTPUT:
% curlZ = Wind stress curl [N/m^3]
%
% DISCLAIMER:
% Albeit this function is designed only for academic purpose, it can be implemented in
% research. Nonetheless, author does not guarantee the accuracy.
%
% REFERENCE:
% A.E. Gill, 1982, “Atmosphere-Ocean Dynamics”, Academy Press, Vol. 30.
% W. G. Large & S. Pond., 1981,“Open Ocean Measurements in Moderate to Strong Winds”,
% J. Physical Oceanography, Vol. 11, pp. 324 - 336.
% K.E. Trenberth, W.G. Large & J.G. Olson, 1990, “The Mean Annual Cycle in Global Ocean
% Wind Stress”, J.Physical Oceanography, Vol. 20, pp. 1742 – 1760.
% NRSC, 2013, "OSCAT Wind stress and Wind stress curl products", Ocean Sciences Group,
% Earth and Climate Science Area, Hyderabad, India.