anchor.diff:  ! 		then compute the plans serially
   ! 		then perform ffts etc in parallel
   ! 		finally destroy plans serially
  +		m=1
           do z=1,nz
               U=sum(domain%u(:realnx-1,z,:))/((realnx-1)*realny)
               V=sum(domain%v(:,z,:realny-1))/(realnx*(realny-1))
               sig  = U*k+V*l
  -            where(sig==0.0) sig=1e-10
  +            where(sig==0.0) sig=1e-15
               denom = sig**2!-f**2
   			
   ! 			where(denom.eq.0) denom=1e-20
  
positive.code: 		logical, intent(in), optional :: reverse
          complex,parameter :: j= (0,1)
          real::gain,offset !used in setting up k and l arrays
          integer::nx,ny,nz,i,midpoint,z,realnx,realny
          real::U,V
  		real,parameter::pi=3.1415927
          type(C_PTR) :: plan
  
negative.code: 	        allocate(vhat(nx,ny))
  	        allocate(v_hat(nx,ny))
  	        allocate(m(nx,ny))
  	        allocate(msq(nx,ny))
  	        allocate(ineta(nx,ny))
          
  	!         # % Compute 2D k and l wavenumber fields (could be stored in options or domain or something)
  