A system and method for sensing magnetic anomalies uses a gradiometer having at least one pair of triaxial magnetometer-accelerometer (TMA) sensors. Each TMA sensor has X, Y, Z magnetic sensing axes and X, Y, Z acceleration sensing axes that are parallel to one another and to the X, Y, Z magnetic sensing and acceleration axes of all other TMA sensors. Each TMA sensor outputs components (Bx, By, Bz) of local magnetic fields and components (Ax, Ay, Az) of local gravitational acceleration fields. The components (Bx, By, Bz) and (Ax, Ay, Az) output from each TMA sensor are processed to generate motion-compensated components (Bcx, Bcy, Bcz) of local magnetic fields. A difference is generated between the motion-compensated components (Bcx, Bcy, Bcz) for each pair of TMA sensors thereby generating differential vector field components (&Dgr;Bx, &Dgr;By,&Dgr;Bz). For improved accuracy, the differential vector field components (&Dgr;Bx, &Dgr;By, &Dgr;Bz) are adjusted using the local gravitational acceleration field components and motion-compensated local magnetic field components in order to compensate for gradiometer motion. Gradient components are generated using the differential vector field components (&Dgr;Bx, &Dgr;By, &Dgr;Bz). In general, for each of the magnetic sensing axes, the gradient components Gij are defined by (&Dgr;Bx/&Dgr;j, &Dgr;By/&Dgr;j, &Dgr;Bz/&Dgr;j), wherein &Dgr;j is a distance between a pair of TMA sensors relative to a j-th one of the X, Y, Z magnetic sensing axes. A scalar-quantity gradient contraction defined as C 2 = ∑ i , j ⁢   ⁢ ( G ij ) 2 is generated for each pair of TMA sensors. The gradient contraction C2 is a robust, rotationally-invariant quantity that changes monotonically with proximity to a magnetic target.