anchor.diff:                       !--------------------------
                        if (p%NBodyMod == 2) then
                           !> The phase shift due to an (x,y) offset for second order difference frequencies is of the form
  -                        !! \f$  exp[-\imath ( k(\omega_1) ( X cos(\Beta(w_1)) + Y sin(\beta(w_1)) )
  -                        !!                  1 k(\omega_2) ( X cos(\Beta(w_2)) + Y sin(\beta(w_2)) ) ) ]\f$.
  +                        !! \f$  exp[-\imath ( k(\omega_1) ( X cos(\beta(w_1)) + Y sin(\beta(w_1)) )
  +                        !!                  1 k(\omega_2) ( X cos(\beta(w_2)) + Y sin(\beta(w_2)) ) ) ]\f$.
                           !! For the first term, \f$ \omega_1 = \omega_2 \$f.
                           !  NOTE: the phase shift applies to the aWaveElevC of the incoming wave.  Including it here instead
                           !        of above is mathematically equivalent, but only because each frequency has only one wave
  
positive.code:                         !--------------------------
                          if (p%NBodyMod == 2) then
                             !> The phase shift due to an (x,y) offset for second order difference frequencies is of the form
                             !! \f$  exp[-\imath ( k(\omega_1) ( X cos(\Beta(w_1)) + Y sin(\beta(w_1)) )
                             !!                  - k(\omega_2) ( X cos(\Beta(w_2)) + Y sin(\beta(w_2)) ) ) ]\f$
                             !  NOTE: the phase shift applies to the aWaveElevC of the incoming wave.  Including it here instead
                             !        of above is mathematically equivalent, but only because each frequency has only one wave
                             !        direction associated with it through the equal energy approach used in multidirectional waves.
  
negative.code:                         !--------------------------
                          if (p%NBodyMod == 2) then
                             !> The phase shift due to an (x,y) offset for second order difference frequencies is of the form
                             !! \f$  exp[-\imath ( k(\omega_1) ( X cos(\beta(w_1)) + Y sin(\beta(w_1)) )
                             !!                  - k(\omega_2) ( X cos(\beta(w_2)) + Y sin(\beta(w_2)) ) ) ]\f$
                             !  NOTE: the phase shift applies to the aWaveElevC of the incoming wave.  Including it here instead
                             !        of above is mathematically equivalent, but only because each frequency has only one wave
                             !        direction associated with it through the equal energy approach used in multidirectional waves.
  