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)) )
  -                           !!                  - 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)) )
  +                           !!                  - 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.
  
positive.code:             if (p%NBodyMod == 2) then
  
                 !> The phase shift due to an (x,y) offset is of the form
                 !! \f$  exp(-\imath k(\omega) ( X cos(\Beta(w)) + Y sin(\beta(w)) )) \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.
  