The derivation of gradient maps in geophysics, particularly in the field of self-potential has the potential to improve our understanding on the source of a signal. Self-potential/elevation gradient maps are beneficial in significantly reducing the topographic effect. Manual calculation of the gradient for large data sets in two-dimensions is time consuming and highly dependent on the direction of the calculation. Automation of the calculation process has the potential to overcome the time and directional dependency problems. The derivation of gradient maps in the multi-direction improves the result and array based operators can perform the automatic calculations rapidly. Four different gradient calculation methods based on a new automatic array oriented procedure (swirl procedure) are discussed and tested with artificial and field data sets. These four methods can be simply defined by the number of data contributing to the calculation (full-swirl or limited-swirl procedures) and the mathematical operator (maximum value or mean value) used in the calculation. The mean value operator using the full-swirl procedure gave the most reliable result in terms of gradient range and accuracy. The swirl procedure can effectively perform the self-potential/elevation gradient calculations and it has a potential use in various applications.