module RotCig_mod

  use rotcig_types
  use fft_types

  use fft_data_mod

  implicit none

  integer, external ::  sreed,srite,sseek_block

  logical linear_interp
  integer, private :: n_z,n_mx,n_z_pos,n_mx_pos
  integer, private :: n_z_pad,n_mx_pad,n_z_pad_pos,n_mx_pad_pos
  integer, private :: n_hx,n_hz,n_hdip,n_hdip_pad,n_hdip_pad_pos,n_hrot,n_hpar
  integer, private :: n_bytes_panel,n_bytes_cube
  integer, private :: n_bytes_nz_cmpx,n_bytes_hdip_cmpx,n_bytes_hpar_cmpx
  integer, private :: n_bytes_gamma_cmpx
  integer, private :: n_threads,n_gamma
  real, private :: d_z,d_mx,d_hx,d_hz,d_hdip,d_hrot,d_hpar,d_gamma
  real, private :: z_0,mx_0,hx_0,hz_0,hdip_0,hrot_0,hpar_0,gamma_0
  real, private :: kz_0,kmx_0,khdip_0,d_kz,d_kmx,d_khdip
  real, private :: pi,sin_rot,cos_rot
  real, private :: shift_z,shift_mx
  real*8, private :: sign_tan


  real, allocatable, dimension (:), private :: k_mx,k_z,k_hdip
  real, allocatable, dimension (:,:), private :: panel,panel_rot
  real, allocatable, dimension (:,:,:), private :: cube
  complex, allocatable, dimension (:,:), private :: temp_z,temp_x
  complex, allocatable, dimension (:,:), private :: k_panel_out,k_panel_inp
  complex, allocatable, dimension (:,:), private :: k_panel_dip,k_panel_par
  complex, allocatable, dimension (:,:), private :: k_panel_gamma
  complex, allocatable, dimension (:,:), private :: k_row,k_temp
  complex, allocatable, dimension (:,:,:), private :: k_cube
!  complex, allocatable, dimension (:,:,:), private :: k_cube_inp,k_cube_out
  complex, allocatable, dimension (:,:,:), private :: temp_x_hdip
!  type (fft_type), private :: fft_nz_id,fft_nmx_id,fft_nhdip_id

  type (fft_type), dimension(:), pointer, private :: fft_nz_id,fft_nmx_id,fft_nhdip_id
  type (fft_type), dimension(:), pointer, private :: fft_nmx_nhdip_id

contains
!------------------------------------------------------------------
  function RotCig_init(rotcig_par)  result(status)
    type (rotcig_param_type) rotcig_par
    integer status

    integer i_k_mx,i_k_z,i_k_hdip,ierr,i_z,i_thread
    real kz,kmx,khdip

    status=0

    pi=acos(-1.)

    linear_interp=rotcig_par%linear_interp

    shift_z=rotcig_par%rot_center_z- rotcig_par%z_0
    shift_mx=rotcig_par%rot_center_mx- rotcig_par%mx_0


    n_z=rotcig_par%n_z
    n_z_pad=rotcig_par%n_z_pad
    n_z_pos=rotcig_par%n_z/2
    n_z_pad_pos=rotcig_par%n_z_pad/2
    if(mod(n_z_pad,2) /=0) then
      write(0,*)'n_z_pad=',n_z_pad
      call seperr ("n_z_pad must be even")
    end if


    n_mx=rotcig_par%n_mx
    n_mx_pad=rotcig_par%n_mx_pad
    n_mx_pos=rotcig_par%n_mx/2
    n_mx_pad_pos=rotcig_par%n_mx_pad/2
    if(mod(n_mx_pad,2) /=0) then
      write(0,*)'n_mx_pad=',n_mx_pad
      call seperr ("n_mx_pad must be even")
    end if

    n_hx=rotcig_par%n_hx
    n_hz=rotcig_par%n_hz
    n_hrot=rotcig_par%n_hrot
    n_hdip=rotcig_par%n_hdip
    n_hdip_pad=rotcig_par%n_hdip_pad
    n_hdip_pad_pos=(n_hdip_pad/2)
    n_hpar=rotcig_par%n_hpar

    d_z=rotcig_par%d_z
    d_mx=rotcig_par%d_mx
    d_hx=rotcig_par%d_hx
    d_hz=rotcig_par%d_hz
    d_hrot=rotcig_par%d_hrot
    d_hdip=rotcig_par%d_hdip
    d_hpar=rotcig_par%d_hpar

    z_0=rotcig_par%z_0
    mx_0=rotcig_par%mx_0
    hx_0=rotcig_par%hx_0
    hz_0=rotcig_par%hz_0
    hrot_0=rotcig_par%hrot_0
    hdip_0=rotcig_par%hdip_0
    hpar_0=rotcig_par%hpar_0

    n_gamma=rotcig_par%n_ang
    gamma_0=(pi*rotcig_par%ang_0/180.)
    d_gamma=(pi*rotcig_par%d_ang/180.)

    n_bytes_panel=n_z*n_mx*4
    n_bytes_cube=n_z*n_mx*n_hdip*4
    n_bytes_nz_cmpx=(n_z_pad_pos+1)*8
    n_bytes_hdip_cmpx=(n_z_pad_pos+1)*n_hdip*8
    n_bytes_hpar_cmpx=(n_z_pad_pos+1)*n_hpar*8
    n_bytes_gamma_cmpx=(n_z_pad_pos+1)*n_gamma*8

    kz=0.
    kmx_0=0.
    d_kz=2.*pi/(d_z*FLOAT(n_z_pad))
    d_kmx=2.*pi/(d_mx*FLOAT(n_mx_pad))

    n_threads=rotcig_par%n_threads
    if(rotcig_par%create_off_ddcig) then
	      if(mod(n_hdip,n_threads) /=0) then
	        write(0,*)'n_threads=',n_threads
	        write(0,*)'n_hdip=',n_hdip
	        call seperr ("n_hdip be multiple of n_threads")
	      end if
	      if(mod(n_hpar,n_threads) /=0) then
	        write(0,*)'n_threads=',n_threads
	        write(0,*)'n_hpar=',n_hpar
	        call seperr ("n_hpar be multiple of n_threads")
	      end if
	      if((rotcig_par%create_ang_ddcig) .and. (rotcig_par%create_ang_cig)) then
		        if(mod(n_z_pad_pos,n_threads) /=0) then
		          write(0,*)'n_threads=',n_threads
		          write(0,*)'n_z_pad_pos=',n_z_pad_pos
		          call seperr ("n_z_pad_pos be multiple of n_threads")
		        end if
		!        if(mod(n_gamma,n_threads) /=0) then
		!          write(0,*)'n_threads=',n_threads
		!          write(0,*)'n_gamma=',n_gamma
		!          call seperr ("n_gamma be multiple of n_threads")
		!        end if
	    else
		        if(mod(n_hx,n_threads) /=0) then
		          write(0,*)'n_threads=',n_threads
		          write(0,*)'n_hx=',n_hrot
		          call seperr ("n_hx be multiple of n_threads")
		        end if
		        if(mod(n_hz,n_threads) /=0) then
		          write(0,*)'n_threads=',n_threads
		          write(0,*)'n_hz=',n_hrot
		          call seperr ("n_hz be multiple of n_threads")
		        end if
		        if(mod(n_hrot,n_threads) /=0) then
		          write(0,*)'n_threads=',n_threads
		          write(0,*)'n_hrot=',n_hrot
		          call seperr ("n_hrot be multiple of n_threads")
		        end if
	    end if
    end if
	
    allocate (k_z(n_z_pad),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of k_z ")
    allocate (k_mx(n_mx_pad),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of k_mx ")
	
    allocate (panel(n_z,n_mx),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of panel ")
    allocate (panel_rot(n_z,n_mx),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of panel_rot ")

    allocate (temp_z(n_z_pad,n_threads),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of temp_z ")

    allocate (temp_x(n_mx_pad,n_threads),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of temp_x ")

    allocate (k_panel_out(n_z_pad_pos+1,n_mx_pad),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of k_panel_out ")

    if(rotcig_par%create_off_ddcig) then
	      allocate (k_panel_dip(n_z_pad_pos+1,n_hdip),stat=ierr)  
	      if(ierr .ne. 0) call seperr("Problems with allocation of k_panel_dip ")
	      allocate (k_panel_par(n_z_pad_pos+1,n_hpar),stat=ierr)  
	      if(ierr .ne. 0) call seperr("Problems with allocation of k_panel_par ")
    else
	      if(rotcig_par%create_cig_tilt) then
		        if(rotcig_par%positive_rot) then
		          sign_tan=1
		        else
		          sign_tan=-1
		        end if
			write(0,*)'sign_tan=',sign_tan
	
		        allocate (k_hdip(n_hdip_pad),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_hdip ")
		        allocate (temp_x_hdip(n_mx_pad,n_hdip_pad,n_threads),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of temp_x_hdip ")
		        allocate (cube(n_z,n_mx,n_hdip),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of cube ")
		        allocate (k_cube(n_z_pad_pos+1,n_mx_pad,n_hdip_pad),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_cube ")
		        allocate (k_row(n_hdip_pad,n_threads),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_row ")
		        allocate (k_temp(n_hdip_pad,n_threads),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_temp ")
		!        allocate (k_cube_out(n_hdip_pad,n_z_pad_pos+1,n_mx_pad),stat=ierr)  
		!        if(ierr .ne. 0) call seperr("Problems with allocation of k_cube_out ")
	
		        khdip_0=0.
		        d_khdip=2.*pi/(d_hdip*FLOAT(n_hdip_pad))
	
		        do i_k_hdip=1,n_hdip_pad
		           khdip  = i_k_hdip - 1.
		           if(khdip <= float(n_hdip_pad/2)) then
		              k_hdip(i_k_hdip) =  khdip_0+d_khdip* khdip
		           else
		              k_hdip(i_k_hdip) = -khdip_0+d_khdip*(khdip-n_hdip_pad)
		           end if
		        end do
	
		        call snap('K_hdip.H',n_hdip_pad,1,k_hdip)
	
! remap k_hdip to facilitate interpolation
		        k_row(:,1)=k_hdip(:)
		        do i_k_hdip=1,n_hdip_pad_pos+1
		          k_hdip(i_k_hdip+n_hdip_pad_pos-1)=k_row(i_k_hdip,1)
		        end do
		        do i_k_hdip=n_hdip_pad_pos+2,n_hdip_pad
		          k_hdip(i_k_hdip-n_hdip_pad_pos-1)=k_row(i_k_hdip,1)
		        end do
	
		        call snap('K_hdip.H',n_hdip_pad,1,k_hdip)

	      else
		      if((rotcig_par%create_ang_ddcig) .or. (rotcig_par%create_ang_cig)) then

		        allocate (k_panel_dip(n_z_pad_pos+1,n_hdip),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_panel_dip ")

		        allocate (k_row(n_hdip_pad,n_threads),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_row ")
		        allocate (k_temp(n_hdip_pad,n_threads),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_temp ")
		        allocate (k_panel_gamma(n_z_pad_pos+1,n_gamma),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_panel_gamma ")


		      else
		        allocate (k_panel_inp(n_z_pad,n_mx_pad),stat=ierr)  
		        if(ierr .ne. 0) call seperr("Problems with allocation of k_panel_inp ")
		      end if

	      end if
    end if

    do i_k_z=1,n_z_pad
       kz  = i_k_z - 1.
       if(kz <= float(n_z_pad/2)) then
          k_z(i_k_z) =  kz_0+d_kz* kz
       else
          k_z(i_k_z) = -kz_0+d_kz*(kz-n_z_pad)
       end if
    end do

    do i_k_mx=1,n_mx_pad
       kmx  = i_k_mx - 1.
       if(kmx <= float(n_mx_pad/2)) then
          k_mx(i_k_mx) =  kmx_0+d_kmx* kmx
       else
          k_mx(i_k_mx) = -kmx_0+d_kmx*(kmx-n_mx_pad)
       end if
    end do

    allocate (fft_nz_id(n_threads),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of fft_nz_id")
    if(rotcig_par%debug2) write(0,*)'Allocate fft_nz_id with shape=',shape(fft_nz_id)

    do i_thread=1,n_threads
      call init_fft_1d(fft_nz_id(i_thread),n_z_pad, &
              temp_z(:,i_thread),temp_z(:,i_thread),fast_init=.false.)  
      if(rotcig_par%debug1) write(0,*)'scale (n_z_pad)=',fft_nz_id(i_thread)%scale
    end do
!    status=fft_nz_id%ierr

    allocate (fft_nmx_id(n_threads),stat=ierr)  
    if(ierr .ne. 0) call seperr("Problems with allocation of fft_nmx_id")
    if(rotcig_par%debug2) write(0,*)'Allocate fft_nmx_id with shape=',shape(fft_nmx_id)

    do i_thread=1,n_threads
      call init_fft_1d(fft_nmx_id(i_thread),n_mx_pad, &
              temp_x(:,i_thread),temp_x(:,i_thread),fast_init=.false.)  
      if(rotcig_par%debug1) write(0,*)'scale (n_mx_pad)=',fft_nmx_id(i_thread)%scale
    end do
!    status=fft_nmx_id%ierr


    if((rotcig_par%create_ang_ddcig) .or. (rotcig_par%create_ang_cig)) then
	      allocate (fft_nhdip_id(n_threads),stat=ierr)  
	      if(ierr .ne. 0) call seperr("Problems with allocation of fft_nhdip_id")
	      if(rotcig_par%debug2) write(0,*)'Allocate fft_nhdip_id with shape=',shape(fft_nhdip_id)

	      do i_thread=1,n_threads
	        call init_fft_1d(fft_nhdip_id(i_thread),n_hdip_pad, &
	                k_row(:,i_thread),k_row(:,i_thread),fast_init=.false.)  
	        if(rotcig_par%debug1) write(0,*)'scale (n_z_pad)=',fft_nhdip_id(i_thread)%scale
	      end do
	    !    status=fft_nhdip_id%ierr
    end if

    if(rotcig_par%create_cig_tilt) then
	      allocate (fft_nmx_nhdip_id(n_threads),stat=ierr)  
	      if(ierr .ne. 0) call seperr("Problems with allocation of fft_nmx_nhdip_id")
	      if(rotcig_par%debug2) write(0,*)'Allocate fft_nmx_nhdip_id with shape=',shape(fft_nmx_nhdip_id)

	      do i_thread=1,n_threads
	        call init_fft_2d(fft_nmx_nhdip_id(i_thread),n_mx_pad,n_hdip_pad, &
	                temp_x_hdip(:,:,i_thread),temp_x_hdip(:,:,i_thread),fast_init=.false.)  
	        if(rotcig_par%debug1) write(0,*)'scale (n_mx_n_hdip)=',fft_nmx_nhdip_id(i_thread)%scale
	      end do
	    !    status=fft_nhdip_id%ierr
    end if

    sin_rot=sin(pi*rotcig_par%out_dip_center/180.)
    cos_rot=cos(pi*rotcig_par%out_dip_center/180.)

    if(rotcig_par%debug1) then
	write(0,*)'ierr (nmx)=',fft_nmx_id%ierr
	write(0,*)'status=',status
	write(0,*)'n_hz=',n_hz
	write(0,*)'n_hx=',n_hx
	write(0,*)'n_hdip=',n_hdip
	write(0,*)'n_hpar=',n_hpar
	write(0,*)'n_hrot=',n_hrot
	write(0,*)'n_gamma=',n_gamma
	write(0,*)'hdip_0=',hdip_0
	write(0,*)'hpar_0=',hpar_0
	write(0,*)'hrot_0=',hrot_0
	write(0,*)'gamma_0=',gamma_0
	write(0,*)'d_hdip=',d_hdip
	write(0,*)'d_hpar=',d_hpar
	write(0,*)'d_hrot=',d_hrot
	write(0,*)'d_gamma=',d_gamma
	write(0,*)'n_z=',n_z
	write(0,*)'n_mx=',n_mx
	write(0,*)'n_z_pad=',n_z_pad
	write(0,*)'n_mx_pad=',n_mx_pad
	write(0,*)'n_hdip_pad=',n_hdip_pad
	write(0,*)'n_z_pos=',n_z_pos
	write(0,*)'n_mx_pos=',n_mx_pos
	write(0,*)'n_z_pad_pos=',n_z_pad_pos
	write(0,*)'n_mx_pad_pos=',n_mx_pad_pos
	write(0,*)'n_hdip_pad_pos=',n_hdip_pad_pos
	write(0,*)'d_kz=',d_kz
	write(0,*)'d_kmx=',d_kmx
	write(0,*)'d_khdip=',d_khdip
	write(0,*)'n_bytes_panel=',n_bytes_panel
	write(0,*)'n_bytes_nz_cmpx=',n_bytes_nz_cmpx
	write(0,*)'n_bytes_hdip_cmpx=',n_bytes_hdip_cmpx
	write(0,*)'n_bytes_hpar_cmpx=',n_bytes_hpar_cmpx
	write(0,*)'n_bytes_gamma_cmpx=',n_bytes_gamma_cmpx
	write(0,*)'sin_rot=',sin_rot
	write(0,*)'cos_rot=',cos_rot
	write(0,*)'shift_z=',shift_z
	write(0,*)'shift_mx=',shift_mx
     end if
	
  end function RotCig_init
!--------------------------------------------------------------------------
  function RotCig_oper(rotcig_par) result(status)

    type (rotcig_param_type) rotcig_par
    integer status

    logical hz_rot_inside,hx_rot_inside,hdip_inside,inside_k_hdip
    integer ret_sseek
    integer i_hrot,i_hz_rot,i_hx_rot,i_hpar,i_hdip,i_mx,i_z,i_thread
    integer i_k_hdip,i_gamma
    real h_rot,hz_rot,hx_rot,hz_min,hz_max,hx_min,hx_max,kz,kmx,kz_out
    real w_hz,w_hx,w_hz_rot,w_hx_rot,sin_rot_map,cos_rot_map
    real hdip_min,hdip_max,h_par,h_dip
    real tan_dip,dip,diff_dip,sin_diff,cos_diff,w_hdip,z_shift
    real tan_gamma,k_hdip_inp,inv_d_khdip,inv_d_z

    if(rotcig_par%create_off_ddcig) then
    !{

      if(rotcig_par%debug2) write(0,*)'Beginning transforming cig_hdip'
      status= RotCig_transf( &
           .true.,rotcig_par%tag_cig_hdip,rotcig_par%tag_cig_hdip_transf,n_hdip,rotcig_par)
      if(rotcig_par%debug2) write(0,*)'Ending transforming cig_hdip'

      if(rotcig_par%debug2) write(0,*)'Beginning remapping offsets'


      ret_sseek=sseek_block( &
            trim(rotcig_par%tag_cig_hdip_transf),0,n_bytes_hdip_cmpx,0)
      if (ret_sseek .ne. 0) then
      !{
        write(0,*)'ret_sseek=',ret_sseek
        write(0,*)'tag_cig_hdip_transf=',trim(rotcig_par%tag_cig_hdip_transf)
        call seperr("Error block seeking cig_hdip_transf in input")
      !}
      end if

      hdip_min=hdip_0
      hdip_max=hdip_0+(n_hdip-1)*d_hdip
      inv_d_z=1./d_z
  
      do i_mx=1,n_mx_pad
      !{

        if(sreed(trim(rotcig_par%tag_cig_hdip_transf),k_panel_dip,n_bytes_hdip_cmpx)/= &
            n_bytes_hdip_cmpx) then
          write(0,*)'tag_cig_hdip_transf=',trim(rotcig_par%tag_cig_hdip_transf)
          call seperr('Problems with reading tag_cig_hdip_transf file')
        end if

        do i_hpar=1,n_hpar
        !{
          h_par=hpar_0+(i_hpar-1)*d_hpar
          if(rotcig_par%debug2) then
            if(mod(i_mx-2,n_mx_pad/4) ==0 ) then
              write(0,*)'i_mx=',i_mx
              write(0,*)'h_par=',h_par
            end if
          end if
          do i_z=1,n_z_pad_pos+1
          !{

            if(i_z /= 1) then
              tan_dip=k_mx(i_mx)/k_z(i_z)
              dip=-atan(tan_dip)
            else
              dip=.5*pi
            end if
            diff_dip=dip-(pi*rotcig_par%inp_dip_center)/180.
            cos_diff=cos(diff_dip)

            if(abs(cos_diff) > 0.) then
!              h_dip=h_par
              h_dip=h_par/cos_diff
              if(rotcig_par%rotate_off_ddcig) then
                h_dip=h_par/(cos_diff**2)
!                h_dip=h_par/(cos_diff)
!                h_dip=h_par
                sin_diff=sin(diff_dip)
!                sin_diff=0.
!                z_shift=h_par*sin_diff*inv_d_z
                z_shift=h_par*sin_diff
!                z_shift=.5*h_par
!                z_shift=-10*d_hpar*sin_diff
              end if
              hdip_inside=((h_dip >= hdip_min) .and. (h_dip < hdip_max))
            else
              hdip_inside=.false.
            end if

            if(rotcig_par%debug2) then
              if(mod(i_mx-1,n_mx_pad/4) ==0 ) then
                if(mod(i_z-1,n_z_pad_pos/4) ==0 ) then
                  write(0,*)'i_z=',i_z
                  write(0,*)'dip=',dip*180/pi
                  write(0,*)'diff_dip=',diff_dip*180/pi
                  write(0,*)'cos_diff=',cos_diff
                  write(0,*)'h_dip=',h_dip
                  write(0,*)'hdip_inside=',hdip_inside
                  if(rotcig_par%rotate_off_ddcig) then
                    write(0,*)'sin_diff=',sin_diff
                    write(0,*)'z_shift=',z_shift
                  end if
                end if
              end if
            end if

            if(hdip_inside) then
            !{

              i_hdip=floor((h_dip-hdip_0)/d_hdip) +1
              w_hdip=     ((h_dip-hdip_0)/d_hdip) +1  - i_hdip

              if(rotcig_par%debug2) then
                if(mod(i_mx-2,n_mx_pad/4) ==0 ) then
                  if(mod(i_z-2,n_z_pad_pos/4) ==0 ) then
                    write(0,*)'i_hdip=',i_hdip
                    write(0,*)'w_hdip=',w_hdip
                  end if
                end if
              end if

              k_panel_par(i_z,i_hpar)= (1.-w_hdip)*k_panel_dip(i_z,i_hdip  ) + &
                                       (   w_hdip)*k_panel_dip(i_z,i_hdip+1)

              if(rotcig_par%rotate_off_ddcig) then
                k_panel_par(i_z,i_hpar)=k_panel_par(i_z,i_hpar)* &
                          cexp(cmplx(0,k_z(i_z)*z_shift))
              end if

            !}
             else
            !{
              k_panel_par(i_z,i_hpar)= cmplx(0.,0.)
             end if
            !}

          !}
          end do
        !}
        end do


        if(srite(trim(rotcig_par%tag_off_ddcig_transf),k_panel_par,n_bytes_hpar_cmpx)/= &
            n_bytes_hpar_cmpx) then
          write(0,*)'tag_off_ddcig_transf=',trim(rotcig_par%tag_off_ddcig_transf)
          call seperr('Problems with writing tag_off_ddcig_transf file')
        end if

      !}
      end do

      if(rotcig_par%debug2) write(0,*)'Ending remapping offsets'

      if(rotcig_par%debug2) write(0,*)'Beginning transforming off_ddcig'
      status= RotCig_transf( &
           .false.,rotcig_par%tag_off_ddcig_transf,rotcig_par%tag_off_ddcig,n_hpar,rotcig_par)
      if(rotcig_par%debug2) write(0,*)'Ending transforming off_ddcig'

    !}
    else
    !{
      if(rotcig_par%create_cig_tilt) then
      !{

      status=RotCig_tilt( &
         rotcig_par%tag_cig_hx,rotcig_par%tag_cig_tilt,rotcig_par)
  
      !}
      else
      !{
      if((rotcig_par%create_ang_ddcig) .or. (rotcig_par%create_ang_cig)) then
      !{

        inv_d_khdip=1./d_khdip

        if(rotcig_par%debug2) write(0,*)'Beginning transforming cig_hdip'
        status= RotCig_transf( &
             .true.,rotcig_par%tag_cig_hdip,rotcig_par%tag_cig_hdip_transf,n_hdip,rotcig_par)
        if(rotcig_par%debug2) write(0,*)'Ending transforming cig_hdip'


        ret_sseek=sseek_block( &
              trim(rotcig_par%tag_ang_ddcig_transf),0,n_bytes_gamma_cmpx,0)
        if (ret_sseek .ne. 0) then
        !{
          write(0,*)'ret_sseek=',ret_sseek
          write(0,*)'tag_ang_ddcig_transf=',trim(rotcig_par%tag_ang_ddcig_transf)
          call seperr("Error block seeking ang_ddcig_transf in input")
        !}
        end if


        ret_sseek=sseek_block( &
              trim(rotcig_par%tag_cig_hdip_transf),0,n_bytes_hdip_cmpx,0)
        if (ret_sseek .ne. 0) then
        !{
          write(0,*)'ret_sseek=',ret_sseek
          write(0,*)'tag_cig_hdip_transf=',trim(rotcig_par%tag_cig_hdip_transf)
          call seperr("Error block seeking cig_hdip_transf in input")
        !}
        end if

        do i_mx=1,n_mx_pad
        !{
          kmx=k_mx(i_mx)
  
          if(sreed(trim(rotcig_par%tag_cig_hdip_transf),k_panel_dip,n_bytes_hdip_cmpx)/= &
              n_bytes_hdip_cmpx) then
            write(0,*)'tag_cig_hdip_transf=',trim(rotcig_par%tag_cig_hdip_transf)
            call seperr('Problems with reading tag_cig_hdip_transf file')
          end if


! BB potential bug with n_thread >1 
          do i_thread=1,n_threads
          !{
          do i_z=1,n_z_pad_pos+1,n_threads
          !{
            kz=k_z(i_z)

            if(rotcig_par%create_ang_ddcig) then
              if(i_z /= 1) then
                tan_dip=kmx/kz
                dip=-atan(tan_dip)
              else
                dip=.5*pi
              end if
              diff_dip=dip-(pi*rotcig_par%inp_dip_center)/180.
              cos_diff=cos(diff_dip)
              sin_diff=sin(diff_dip)

              kz_out=  kz*cos_diff - kmx*sin_diff
            else
              cos_diff=1
              kz_out=  kz
            end if

          if(rotcig_par%debug2) then
            if(mod(i_mx-2,n_mx_pad/4) ==0 ) then
              if(mod(i_z-2,n_z_pad_pos/4) ==0 ) then
                  write(0,*)'tan_dip=',tan_dip
                  write(0,*)'dip=',dip
                  write(0,*)'cos_diff=',cos_diff
                  write(0,*)'kmx=',kmx
                  write(0,*)'kz=',kz
                  write(0,*)'kz_out=',kz_out
                end if
             end if
          end if


            do i_hdip=1,n_hdip
            !{
              k_row(i_hdip,i_thread)=k_panel_dip(i_z,i_hdip)
            !}
            end do
            do i_hdip=n_hdip+1,n_hdip_pad
            !{
              k_row(i_hdip,i_thread)=cmplx(0.0)
            !}
            end do
  
!            call fft_1d_data(fft_nhdip_id(i_thread),FFT_FORW,n_hdip_pad,k_row(:,i_thread))
            call fft_1d_data(fft_nhdip_id(i_thread),FFT_BACK,n_hdip_pad,k_row(:,i_thread))
            do i_hdip=1,n_hdip_pad_pos+1
            !{
              k_temp(i_hdip+n_hdip_pad_pos-1,i_thread)=k_row(i_hdip,i_thread)
            !}
            end do
            do i_hdip=n_hdip_pad_pos+2,n_hdip_pad
            !{
              k_temp(i_hdip-n_hdip_pad_pos-1,i_thread)=k_row(i_hdip,i_thread)
            !}
            end do

            do i_gamma=1,n_gamma
            !{
  
              tan_gamma=(gamma_0 + (i_gamma-1)*d_gamma)

              k_hdip_inp=kz_out*tan_gamma*cos_diff

              i_k_hdip=floor(k_hdip_inp*inv_d_khdip+float(n_hdip_pad_pos))
              inside_k_hdip=((i_k_hdip >= 1) .and.(i_k_hdip < n_hdip_pad))

              w_hdip=(k_hdip_inp*inv_d_khdip)+n_hdip_pad_pos - i_k_hdip

!write(0,*)'#9'

              if(inside_k_hdip) then
                k_panel_gamma(i_z,i_gamma)= &
                         (1.-w_hdip)*k_temp(i_k_hdip  ,i_thread) + &
                         (   w_hdip)*k_temp(i_k_hdip+1,i_thread)
!                k_panel_gamma(i_z,i_gamma)= cmplx(gamma,tan_gamma)
              end if

            !}
            end do
          !}
          end do
          !}
          end do
  
          if(srite(trim(rotcig_par%tag_ang_ddcig_transf),k_panel_gamma,n_bytes_gamma_cmpx)/= &
              n_bytes_gamma_cmpx) then
            write(0,*)'tag_ang_ddcig_transf=',trim(rotcig_par%tag_ang_ddcig_transf)
            call seperr('Problems with writing tag_ang_ddcig_transf file')
          end if


        !}
        end do
  
        if(rotcig_par%debug2) write(0,*)'Beginning transforming back ang_ddcig'

        status= RotCig_transf( &
             .false.,rotcig_par%tag_ang_ddcig_transf,rotcig_par%tag_ang_ddcig,n_gamma,rotcig_par)
        if(rotcig_par%debug2) write(0,*)'Ending transforming back ang_ddcig'


      !}
      else
      !{
        if(rotcig_par%cig_hz) then
        !{
          if(rotcig_par%debug2) write(0,*)'Beginning rotating cig_hz'
          status= RotCig_rot_weight( &
               rotcig_par%tag_cig_hz,rotcig_par%tag_cig_hz_rot_weight,n_hz,rotcig_par)
          if(rotcig_par%debug2) write(0,*)'Ending rotating cig_hz'
    
        !}
        end if
        if(rotcig_par%cig_hx) then
        !{
          if(rotcig_par%debug2) write(0,*)'Beginning rotating cig_hx'
          status= RotCig_rot_weight( &
               rotcig_par%tag_cig_hx,rotcig_par%tag_cig_hx_rot_weight,n_hx,rotcig_par)
          if(rotcig_par%debug2) write(0,*)'Ending rotating cig_hx'
    
        !}
        end if
  
        if(rotcig_par%merge_hx_hz) then
        !{
  
          ret_sseek=sseek_block( &
                trim(rotcig_par%tag_cig_hz_rot_weight),0,n_bytes_panel,0)
          if (ret_sseek .ne. 0) then
          !{
            write(0,*)'ret_sseek=',ret_sseek
            write(0,*)'tag_cig_hz_rot_weight=',trim(rotcig_par%tag_cig_hz_rot_weight)
            call seperr("Error block seeking cig_hz_rot_weight in input")
          !}
          end if
        
          ret_sseek=sseek_block( &
                trim(rotcig_par%tag_cig_hx_rot_weight),0,n_bytes_panel,0)
          if (ret_sseek .ne. 0) then
          !{
            write(0,*)'ret_sseek=',ret_sseek
            write(0,*)'tag_cig_hx_rot_weight=',trim(rotcig_par%tag_cig_hx_rot_weight)
            call seperr("Error block seeking cig_hx_rot_weight in input")
          !}
          end if
    
    
          do i_hrot=1,n_hrot
          !{
             if(sreed(trim(rotcig_par%tag_cig_hz_rot_weight),panel,n_bytes_panel)/= &
                 n_bytes_panel) then
               write(0,*)'tag_cig_hz_rot_weight=',trim(rotcig_par%tag_cig_hz_rot_weight)
               call seperr('Problems with reading cig_hz_rot_weight file')
             end if
      
             panel_rot=panel
      
             if(sreed(trim(rotcig_par%tag_cig_hx_rot_weight),panel,n_bytes_panel)/= &
                 n_bytes_panel) then
               write(0,*)'tag_cig_hx_rot_weight=',trim(rotcig_par%tag_cig_hx_rot_weight)
               call seperr('Problems with reading cig_hx_rot_weight file')
             end if
    
             panel_rot=panel_rot + panel
    
            if(srite(trim(rotcig_par%tag_cig_hrot_weight),panel_rot,n_bytes_panel)/= &
                n_bytes_panel) then
              write(0,*)'tag_cig_hrot_weight=',trim(rotcig_par%tag_cig_hrot_weight)
              call seperr('Problems with writing cig_hrot_weight file')
            end if
      
          !}
          end do 
  
        !}
        end if

      !}
      end if
      !}
      end if

    !}
    end if
    
      status=0
    return	
  end function RotCig_oper
!

  function RotCig_rot_weight(tag_cig,tag_cig_rot_weight,n_h,rotcig_par) result(status)

    character(len=128) :: tag_cig,tag_cig_rot_weight
    integer n_h
    type (rotcig_param_type) rotcig_par
    integer status

    logical inside_kz,inside_kmx,cig_hx
    integer i_z,i_mx,i_h,i_kz_inp,i_kmx_inp,i_thread
    real kz,kmx,kz_inp,kmx_inp,inv_d_kz,inv_d_kmx
    real w_z,w_mx,w_cig,dip,tan_dip,epsi

    inv_d_kz=1./d_kz
    inv_d_kmx=1./d_kmx
    epsi=.0001*d_kz

    if(trim(tag_cig) == trim(rotcig_par%tag_cig_hx)) then
      cig_hx=.true.
    else
      if(trim(tag_cig) == trim(rotcig_par%tag_cig_hz)) then
        cig_hx=.false.
      else
        call seperr('RotCig_rot_weight: Unrecognized tag_cig')
      end if
    end if

write(0,*)'cig_hx=',cig_hx
write(0,*)'n_h=',n_h

    do i_thread=1,n_threads
    !{
    do i_h=1,n_h,n_threads
    !{
      if(sreed(trim(tag_cig),panel,n_bytes_panel)/= &
          n_bytes_panel) then
        write(0,*)'tag_cig=',trim(tag_cig)
        call seperr('Problems with reading tag_cig file')
      end if

! FFT along z
      do i_mx=1,n_mx
      !{
        temp_z(1:n_z,i_thread)=cmplx(panel(1:n_z,i_mx),0.)
        do i_z=n_z+1,n_z_pad
          temp_z(i_z,i_thread)=cmplx(0.,0.)
        end do
!        call fft_1d_data(fft_nz_id(i_thread),FFT_FORW,n_z_pad,temp_z(:,i_thread))
        call fft_1d_data(fft_nz_id(i_thread),FFT_BACK,n_z_pad,temp_z(:,i_thread))

        do i_z=1,n_z_pad_pos+1
          k_panel_inp(i_z,i_mx)=temp_z(i_z,i_thread)
        end do

      !}
      end do 
      do i_mx=n_mx+1,n_mx_pad
      !{
        do i_z=1,n_z_pad_pos+1
          k_panel_inp(i_z,i_mx)=cmplx(0.,0.)
        end do
      !}
      end do 

      do i_z=n_z_pad_pos+2,n_z_pad
      !{
        k_panel_inp(i_z,:)=conjg(k_panel_inp(n_z_pad-i_z+2,:))
      !}
      end do

! FFT along mx
      do i_z=1,n_z_pad
      !{
        temp_x(:,i_thread)=k_panel_inp(i_z,:)
!        call fft_1d_data(fft_nmx_id(i_thread),FFT_FORW,n_mx_pad,temp_x(:,i_thread))
        call fft_1d_data(fft_nmx_id(i_thread),FFT_BACK,n_mx_pad,temp_x(:,i_thread))
        k_panel_inp(i_z,:)=temp_x(:,i_thread)
!shift input
        do i_mx=1,n_mx_pad
        !{
          k_panel_inp(i_z,i_mx)=k_panel_inp(i_z,i_mx) * &
           cexp(cmplx(0,-k_z(i_z)*shift_z-k_mx(i_mx)*shift_mx))
        !}
        end do 
      !}
      end do 


! remap input 
      do i_mx=1,n_mx_pad
      !{
        do i_z=1,n_z_pad_pos+1
          temp_z(i_z+n_z_pad_pos-1,i_thread)=k_panel_inp(i_z,i_mx)
        end do
        do i_z=n_z_pad_pos+2,n_z_pad
          temp_z(i_z-n_z_pad_pos-1,i_thread)=k_panel_inp(i_z,i_mx)
        end do
        k_panel_inp(1:n_z_pad,i_mx)=temp_z(1:n_z_pad,i_thread)
      !}
      end do
      do i_z=1,n_z_pad
      !{
        do i_mx=1,n_mx_pad_pos+1
          temp_x(i_mx+n_mx_pad_pos-1,i_thread)=k_panel_inp(i_z,i_mx)
        end do
        do i_mx=n_mx_pad_pos+2,n_mx_pad
          temp_x(i_mx-n_mx_pad_pos-1,i_thread)=k_panel_inp(i_z,i_mx)
        end do
        k_panel_inp(i_z,1:n_mx_pad)=temp_x(1:n_mx_pad,i_thread)
      !}
      end do

! rotate
      do i_mx=1,n_mx_pad
      !{
        kmx=k_mx(i_mx)
        do i_z=1,n_z_pad_pos+1
        !{
          kz=k_z(i_z)


          kz_inp=  kz*cos_rot + kmx*sin_rot
          kmx_inp=-kz*sin_rot + kmx*cos_rot

          if(abs(kz_inp) >= epsi) then
            tan_dip=kmx_inp/kz_inp
            dip=-atan(tan_dip)
          else
            dip=.5*pi
          end if
          if(rotcig_par%merge_hx_hz) then
            if(cig_hx) then
              w_cig=cos(dip)**2
            else
              w_cig=sin(dip)**2
            end if
          else
            w_cig=1.
          end if

          if(rotcig_par%debug2) then
            if(mod(i_mx-2,n_mx_pad/4) ==0 ) then
              if(mod(i_z-2,n_z_pad_pos/4) ==0 ) then
                  write(0,*)'tan_dip=',tan_dip
                  write(0,*)'dip=',dip
                  write(0,*)'w_cig=',w_cig
                end if
             end if
          end if

          if(linear_interp) then
          !{
            i_kz_inp=floor(kz_inp*inv_d_kz+float(n_z_pad_pos))
            i_kmx_inp=floor(kmx_inp*inv_d_kmx+float(n_mx_pad_pos))

            inside_kz=((i_kz_inp >= 1) .and.(i_kz_inp < n_z_pad))
            inside_kmx=((i_kmx_inp >= 1) .and.(i_kmx_inp < n_mx_pad))

            w_z=(kz_inp*inv_d_kz)+n_z_pad_pos - i_kz_inp
            w_mx=(kmx_inp*inv_d_kmx)+n_mx_pad_pos - i_kmx_inp

            if(inside_kz .and. inside_kmx) then
            !{

              k_panel_out(i_z,i_mx)= w_cig*( &
                (1-w_z)*(1-w_mx)*k_panel_inp(i_kz_inp  ,i_kmx_inp  ) + &
                (1-w_z)*(  w_mx)*k_panel_inp(i_kz_inp  ,i_kmx_inp+1) + &
                (  w_z)*(1-w_mx)*k_panel_inp(i_kz_inp+1,i_kmx_inp  ) + &
                (  w_z)*(  w_mx)*k_panel_inp(i_kz_inp+1,i_kmx_inp+1) &
                                           )

            !}
            else
            !{
              if((i_kz_inp == 0) .and.(i_kmx_inp == 0)) then
              !{
                k_panel_out(i_z,i_mx)= w_cig*( &
                  (  w_z)*(  w_mx)*k_panel_inp(i_kz_inp+1,i_kmx_inp+1  ) &
                                             )
              !}
              else
              !{
                if((i_kz_inp == n_z_pad) .and.(i_kmx_inp == 0)) then
                !{
                  k_panel_out(i_z,i_mx)= w_cig*( &
                    (1-w_z)*(  w_mx)*k_panel_inp(i_kz_inp  ,i_kmx_inp+1)  &
                                               )
                !}
                else
                !{
                  if((i_kz_inp == 0) .and.(i_kmx_inp == n_mx_pad)) then
                  !{
                    k_panel_out(i_z,i_mx)= w_cig*( &
                      (  w_z)*(1-w_mx)*k_panel_inp(i_kz_inp+1,i_kmx_inp  )  &
                                                 )
                  !}
                  else
                  !{
                    if((i_kz_inp == n_z_pad) .and.(i_kmx_inp == n_mx_pad))then
                      k_panel_out(i_z,i_mx)= w_cig*( &
                        (1-w_z)*(1-w_mx)*k_panel_inp(i_kz_inp  ,i_kmx_inp  ) &
                                                   )
                    else
                       k_panel_out(i_z,i_mx)=cmplx(0,0)
                    end if
                  !}
                  end if
                !}
                end if
              !}
              end if
            !}
            end if

          !}
          else
          !{
            i_kz_inp=nint(kz_inp*inv_d_kz+float(n_z_pad_pos))
            i_kmx_inp=nint(kmx_inp*inv_d_kmx+float(n_mx_pad_pos))

            inside_kz=((i_kz_inp >= 1) .and.(i_kz_inp <= n_z_pad))
            inside_kmx=((i_kmx_inp >= 1) .and.(i_kmx_inp <= n_mx_pad))

            if(inside_kz .and. inside_kmx) then
               k_panel_out(i_z,i_mx)=w_cig*k_panel_inp(i_kz_inp,i_kmx_inp) 
            else
              k_panel_out(i_z,i_mx)=cmplx(0,0)
            end if
          !}
          end if

!shift output
            k_panel_out(i_z,i_mx)=k_panel_out(i_z,i_mx) * &
             cexp(cmplx(0,k_z(i_z)*shift_z+k_mx(i_mx)*shift_mx))

        !}
        end do
      !}
      end do 

! INV FFT along mx

      do i_z=1,n_z_pad_pos+1
      !{
        temp_x(1:n_mx_pad,i_thread)=k_panel_out(i_z,1:n_mx_pad)
!        call fft_1d_data(fft_nmx_id(i_thread),FFT_BACK,n_mx_pad,temp_x(:,i_thread))
        call fft_1d_data(fft_nmx_id(i_thread),FFT_FORW,n_mx_pad,temp_x(:,i_thread))
        k_panel_out(i_z,1:n_mx_pad)=temp_x(1:n_mx_pad,i_thread)
      !}
      end do 

! NV FFT along z
      do i_mx=1,n_mx
      !{
        do i_z=1,n_z_pad_pos+1
        !{
          temp_z(i_z,i_thread)=k_panel_out(i_z,i_mx)
        !}
        end do
        do i_z=n_z_pad_pos+2,n_z_pad
        !{
          temp_z(i_z,i_thread)= conjg(k_panel_out(n_z_pad-i_z+2,i_mx))
        !}
        end do

!        call fft_1d_data(fft_nz_id(i_thread),FFT_BACK,n_z_pad,temp_z(:,i_thread))
        call fft_1d_data(fft_nz_id(i_thread),FFT_FORW,n_z_pad,temp_z(:,i_thread))

        panel(1:n_z,i_mx)=real(temp_z(1:n_z,i_thread))

      !}
      end do 


      if(srite(trim(tag_cig_rot_weight),panel,n_bytes_panel)/= &
          n_bytes_panel) then
        write(0,*)'tag_cig_rot_weight=',trim(tag_cig_rot_weight)
        call seperr('Problems with writing tag_cig_rot_weight file')
      end if

!write(0,*)' Done i_h=',i_h

    !}
    end do 
    !}
    end do 


    status=0
    return	

  end function RotCig_rot_weight


  function RotCig_transf(fft_forward,tag_cig_inp,tag_cig_out,n_h_gamma,rotcig_par) &
              result(status)

    logical fft_forward
    character(len=128) :: tag_cig_inp,tag_cig_out
    integer n_h_gamma
    type (rotcig_param_type) rotcig_par
    integer status

    integer i_z,i_mx,i_h_gamma,i_thread
    integer ret_sseek,i_sseek


    do i_thread=1,n_threads
    !{
    do i_h_gamma=1,n_h_gamma,n_threads
    !{


      if(fft_forward) then
      !{
        if(sreed(trim(tag_cig_inp),panel,n_bytes_panel)/= &
            n_bytes_panel) then
          write(0,*)'tag_cig=',trim(tag_cig_inp)
          call seperr('Problems with reading tag_cig_inp file')
        end if
  
  ! FFT along z
        do i_mx=1,n_mx
        !{
          temp_z(1:n_z,i_thread)=cmplx(panel(1:n_z,i_mx),0.)
          do i_z=n_z+1,n_z_pad
            temp_z(i_z,i_thread)=cmplx(0.,0.)
          end do
!          call fft_1d_data(fft_nz_id(i_thread),FFT_FORW,n_z_pad,temp_z(:,i_thread))
          call fft_1d_data(fft_nz_id(i_thread),FFT_BACK,n_z_pad,temp_z(:,i_thread))
  
          do i_z=1,n_z_pad_pos+1
            k_panel_out(i_z,i_mx)=temp_z(i_z,i_thread)
          end do
  
        !}
        end do 
        do i_mx=n_mx+1,n_mx_pad
        !{
          do i_z=1,n_z_pad_pos+1
            k_panel_out(i_z,i_mx)=cmplx(0.,0.)
          end do
        !}
        end do 
  
  ! FFT along mx
        do i_z=1,n_z_pad_pos+1
        !{
          temp_x(:,i_thread)=k_panel_out(i_z,:)
!          call fft_1d_data(fft_nmx_id(i_thread),FFT_FORW,n_mx_pad,temp_x(:,i_thread))
          call fft_1d_data(fft_nmx_id(i_thread),FFT_BACK,n_mx_pad,temp_x(:,i_thread))
          k_panel_out(i_z,:)=temp_x(:,i_thread)
        !}
        end do 

        do i_mx=1,n_mx_pad
        !{

          i_sseek=((i_mx-1)*n_h_gamma + (i_h_gamma-1))
          ret_sseek=sseek_block( &
                trim(tag_cig_out),i_sseek,n_bytes_nz_cmpx,0)
          if (ret_sseek .ne. i_sseek) then
          !{
            write(0,*)'ret_sseek=',ret_sseek
            write(0,*)'i_sseek=',i_sseek
            write(0,*)'tag_cig_out=',trim(tag_cig_out)
            call seperr("Error block seeking cig_out in output")
          !}
          end if
    
          do i_z=1,n_z_pad_pos+1
            temp_z(i_z,i_thread)=k_panel_out(i_z,i_mx)
          end do

          if(srite(trim(tag_cig_out),temp_z,n_bytes_nz_cmpx)/= &
              n_bytes_nz_cmpx) then
            write(0,*)'tag_cig_out=',trim(tag_cig_out)
            call seperr('Problems with writing tag_cig_out file')
          end if

        !}
        end do 

      !}
      else
      !{

        do i_mx=1,n_mx_pad
        !{

          i_sseek=((i_mx-1)*n_h_gamma + (i_h_gamma-1))
          ret_sseek=sseek_block( &
                trim(tag_cig_inp),i_sseek,n_bytes_nz_cmpx,0)
          if (ret_sseek .ne. i_sseek) then
          !{
            write(0,*)'ret_sseek=',ret_sseek
            write(0,*)'i_sseek=',i_sseek
            write(0,*)'tag_cig_out=',trim(tag_cig_inp)
            call seperr("Error block seeking cig_inp in output")
          !}
          end if

          if(sreed(trim(tag_cig_inp),temp_z,n_bytes_nz_cmpx)/= &
              n_bytes_nz_cmpx) then
            write(0,*)'tag_cig_inp=',trim(tag_cig_inp)
            call seperr('Problems with reading tag_cig_inp file')
          end if

          do i_z=1,n_z_pad_pos+1
            k_panel_out(i_z,i_mx)=temp_z(i_z,i_thread)
          end do

        !}
        end do 

! INV FFT along mx

        do i_z=1,n_z_pad_pos+1
        !{
          temp_x(1:n_mx_pad,i_thread)=k_panel_out(i_z,1:n_mx_pad)
!          call fft_1d_data(fft_nmx_id(i_thread),FFT_BACK,n_mx_pad,temp_x(:,i_thread))
          call fft_1d_data(fft_nmx_id(i_thread),FFT_FORW,n_mx_pad,temp_x(:,i_thread))
          k_panel_out(i_z,1:n_mx_pad)=temp_x(1:n_mx_pad,i_thread)
        !}
        end do 

! NV FFT along z
        do i_mx=1,n_mx
        !{
          do i_z=1,n_z_pad_pos+1
          !{
            temp_z(i_z,i_thread)=k_panel_out(i_z,i_mx)
          !}
          end do
          do i_z=n_z_pad_pos+2,n_z_pad
          !{
            temp_z(i_z,i_thread)= conjg(k_panel_out(n_z_pad-i_z+2,i_mx))
          !}
          end do
  
!          call fft_1d_data(fft_nz_id(i_thread),FFT_BACK,n_z_pad,temp_z(:,i_thread))
          call fft_1d_data(fft_nz_id(i_thread),FFT_FORW,n_z_pad,temp_z(:,i_thread))
  
          panel(1:n_z,i_mx)=real(temp_z(1:n_z,i_thread))
  
        !}
        end do 


        if(srite(trim(tag_cig_out),panel,n_bytes_panel)/= &
            n_bytes_panel) then
          write(0,*)'tag_cig_out=',trim(tag_cig_out)
          call seperr('Problems with writing tag_cig_out file')
        end if

      !}
      end if

    !}
    end do 
    !}
    end do 


    status=0
    return	

  end function RotCig_transf

  function RotCig_tilt(tag_cig_inp,tag_cig_out,rotcig_par) &
              result(status)

    character(len=128) :: tag_cig_inp,tag_cig_out
    type (rotcig_param_type) rotcig_par
    integer status

    logical inside_k_hdip
    integer i_z,i_mx,i_hdip,i_thread,i_k_hdip
    integer ret_sseek,i_sseek
    real geol_dip,off_dip,inp_dip,shift
    real*8 tan_geol_dip,tan_off_dip,tan_inp_dip,tan_out_dip
    real inv_d_khdip,k_hdip_inp,k_hdip_out,w_hdip

    i_thread=1
    inv_d_khdip=1./d_khdip

    if(sreed(trim(tag_cig_inp),cube,n_bytes_cube)/= &
        n_bytes_cube) then
      write(0,*)'tag_cig=',trim(tag_cig_inp)
      call seperr('Problems with reading tag_cig_inp file')
    end if
  
! zeroes k_cube
    k_cube=cmplx(0.,0.)
    k_temp=cmplx(0.,0.)
    k_row=cmplx(0.,0.)

  ! FFT along z
    do i_mx=1,n_mx
    !{
      do i_hdip=1,n_hdip
      !{

        temp_z(1:n_z,i_thread)=cmplx(cube(1:n_z,i_mx,i_hdip),0.)
        do i_z=n_z+1,n_z_pad
          temp_z(i_z,i_thread)=cmplx(0.,0.)
        end do

        call fft_1d_data(fft_nz_id(i_thread),FFT_BACK,n_z_pad,temp_z(:,i_thread))

        do i_z=1,n_z_pad_pos+1
          k_cube(i_z,i_mx,i_hdip)=temp_z(i_z,i_thread)
        end do

      !}
      end do 
    !}
    end do

  ! FFT along mx and hdip
    do i_z=1,n_z_pad_pos+1
    !{
      temp_x_hdip(:,:,i_thread)=k_cube(i_z,:,:)
      call fft_2d_data(fft_nmx_nhdip_id(i_thread),FFT_BACK,n_mx_pad,n_hdip_pad,temp_x_hdip(:,:,i_thread))
      k_cube(i_z,:,:)=temp_x_hdip(:,:,i_thread)
    !}
    end do 

    do i_mx=1,n_mx_pad
    !{
      do i_z=1,n_z_pad_pos+1
      !{
        if(i_z /= 1) then
          tan_geol_dip=-k_mx(i_mx)/k_z(i_z)
!          geol_dip=-atan(tan_dip)
!          geol_dip=0.
        else
          tan_geol_dip=sign_tan*sign(huge(tan_geol_dip),-k_mx(i_mx))
        end if

        do i_hdip=1,n_hdip_pad_pos+1
        !{
          k_temp(i_hdip+n_hdip_pad_pos-1,i_thread)=k_cube(i_z,i_mx,i_hdip)
        !}
        end do
        do i_hdip=n_hdip_pad_pos+2,n_hdip_pad
        !{
          k_temp(i_hdip-n_hdip_pad_pos-1,i_thread)=k_cube(i_z,i_mx,i_hdip)
        !}
        end do

        do i_hdip=1,n_hdip_pad
        !{
          shift=k_hdip(i_hdip)*hdip_0
          k_temp(i_hdip,i_thread)=k_temp(i_hdip,i_thread)*cexp(cmplx(0,shift))
        !}
        end do

        if(rotcig_par%push_interp) then
        !{
          do i_hdip=1,n_hdip_pad
          !{
            if(i_z /= 1) then
            !{
              tan_off_dip=-sign_tan*k_hdip(i_hdip)/k_z(i_z)
!              if (abs(1-sign_tan*tan_off_dip*tan_geol_dip)>tiny(tan_off_dip)) &
              if (abs(1-tan_off_dip*tan_geol_dip)>tiny(tan_off_dip)) &
                  then
                tan_out_dip= &
                  (tan_geol_dip + tan_off_dip)/&
                  (1 - tan_off_dip*tan_geol_dip)
!                  (tan_geol_dip + sign_tan*tan_off_dip)/&
!                  (1 - sign_tan*tan_off_dip*tan_geol_dip)
  
!                k_hdip_out=-k_z(i_z)*tan_out_dip
                k_hdip_out=-sign_tan*k_z(i_z)*tan_out_dip
              else
                k_hdip_out=sign_tan*sign(k_hdip(i_hdip),-k_mx(i_mx))
              end if
            !}
            else
            !{
                k_hdip_out=sign_tan*sign(k_hdip(i_hdip),-k_mx(i_mx))
            !}
            end if

            i_k_hdip=floor(k_hdip_out*inv_d_khdip+float(n_hdip_pad_pos))
            inside_k_hdip=((i_k_hdip >= 1) .and.(i_k_hdip < n_hdip_pad))
  
            w_hdip=(k_hdip_out*inv_d_khdip)+n_hdip_pad_pos - i_k_hdip

            if(inside_k_hdip) then
              k_row(i_k_hdip,i_thread)=   (1.-w_hdip)*k_temp(i_hdip,i_thread) 
              k_row(i_k_hdip+1,i_thread)= (   w_hdip)*k_temp(i_hdip,i_thread)
            else

              if(i_k_hdip == n_hdip_pad) then
                k_row(i_k_hdip,i_thread)= (1.-w_hdip)*k_temp(i_hdip,i_thread) 
                k_row(       1,i_thread)= (   w_hdip)*k_temp(i_hdip,i_thread)
              else
                if(i_k_hdip == 0) then
                  k_row(n_hdip_pad,i_thread)=(1.-w_hdip)*k_temp(i_hdip,i_thread) 
                  k_row(         1,i_thread)=(   w_hdip)*k_temp(i_hdip,i_thread)
                end if
!                k_row(i_hdip,i_thread)= k_temp(i_k_hdip  ,i_thread) 
              end if
            end if
          !}
          end do 
        !}
        else
        !{

          do i_hdip=1,n_hdip_pad
          !{
            if(i_z /= 1) then
            !{
              tan_off_dip=-k_hdip(i_hdip)/k_z(i_z)
! Notice opposite sign in adding tangents
              if (abs(1+sign_tan*tan_off_dip*tan_geol_dip)>tiny(tan_off_dip)) &
                   then
                tan_inp_dip= &
                  (tan_geol_dip - sign_tan*tan_off_dip)/&
                  (1+sign_tan*tan_off_dip*tan_geol_dip)
  
                k_hdip_inp=-k_z(i_z)*tan_inp_dip
              else
                k_hdip_inp=sign(k_hdip(i_hdip),-k_mx(i_mx))
              end if
            !}
            else
            !{
                k_hdip_inp=sign(k_hdip(i_hdip),-k_mx(i_mx))
            !}
            end if

!            diff_dip=dip-(pi*rotcig_par%inp_dip_center)/180.
!            cos_diff=cos(diff_dip)
!          shift=k_z(i_z)*(1*d_z) - k_mx(i_mx)*(d_mx*10) + & 
!                k_hdip(i_hdip)*(d_hdip*1.)
!          k_row(i_hdip,i_thread)=k_temp(i_hdip,i_thread)* &
!                          cexp(cmplx(0,shift))


!          k_hdip_inp=tan(atan(k_hdip(i_hdip)))
!          k_hdip_inp=k_hdip(i_hdip)

            i_k_hdip=floor(k_hdip_inp*inv_d_khdip+float(n_hdip_pad_pos))
            inside_k_hdip=((i_k_hdip >= 1) .and.(i_k_hdip < n_hdip_pad))
  
            w_hdip=(k_hdip_inp*inv_d_khdip)+n_hdip_pad_pos - i_k_hdip
  
            if(inside_k_hdip) then
              k_row(i_hdip,i_thread)= &
                           (1.-w_hdip)*k_temp(i_k_hdip  ,i_thread) + &
                           (   w_hdip)*k_temp(i_k_hdip+1,i_thread)
            else
              if(i_k_hdip == n_hdip_pad) then
                k_row(i_hdip,i_thread)= &
                             (1.-w_hdip)*k_temp(i_k_hdip  ,i_thread) + &
                             (   w_hdip)*k_temp(1         ,i_thread)
              else
                if(i_k_hdip == 0) then
                  k_row(i_hdip,i_thread)= &
                             (1.-w_hdip)*k_temp(n_hdip_pad,i_thread) + &
                             (   w_hdip)*k_temp(1         ,i_thread)
                end if
!                k_row(i_hdip,i_thread)= k_temp(i_k_hdip  ,i_thread) 
              end if
            end if
  
          !}
          end do 
        !}
        end if

        do i_hdip=1,n_hdip_pad
        !{
          shift=-k_hdip(i_hdip)*hdip_0
          k_row(i_hdip,i_thread)=k_row(i_hdip,i_thread)*cexp(cmplx(0,shift))
        !}
        end do

        do i_hdip=1,n_hdip_pad_pos+1
        !{
          k_cube(i_z,i_mx,i_hdip)=k_row(i_hdip+n_hdip_pad_pos-1,i_thread)
        !}
        end do
        do i_hdip=n_hdip_pad_pos+2,n_hdip_pad
        !{
          k_cube(i_z,i_mx,i_hdip)=k_row(i_hdip-n_hdip_pad_pos-1,i_thread)
        !}
        end do

      !}
      end do 
    !}
    end do 


! INV FFT along mx and hdip

    do i_z=1,n_z_pad_pos+1
    !{
      temp_x_hdip(:,:,i_thread)=k_cube(i_z,:,:)
      call fft_2d_data(fft_nmx_nhdip_id(i_thread),FFT_FORW,n_mx_pad,n_hdip_pad,temp_x_hdip(:,:,i_thread))

       k_cube(i_z,:,:)=temp_x_hdip(:,:,i_thread)
    !}
    end do 

! NV FFT along z
    do i_mx=1,n_mx
    !{
      do i_hdip=1,n_hdip
      !{
        do i_z=1,n_z_pad_pos+1
        !{
          temp_z(i_z,i_thread)=k_cube(i_z,i_mx,i_hdip)
        !}
        end do
        do i_z=n_z_pad_pos+2,n_z_pad
        !{
          temp_z(i_z,i_thread)= conjg(k_cube(n_z_pad-i_z+2,i_mx,i_hdip))
        !}
        end do
!        temp_z(1,i_thread)=cmplx(0.,0.)

        call fft_1d_data(fft_nz_id(i_thread),FFT_FORW,n_z_pad,temp_z(:,i_thread))
        cube(1:n_z,i_mx,i_hdip)=real(temp_z(1:n_z,i_thread))

      !}
      end do 
    !}
    end do 

    if(srite(trim(tag_cig_out),cube,n_bytes_cube)/= &
        n_bytes_cube) then
      write(0,*)'tag_cig_out=',trim(tag_cig_out)
      call seperr('Problems with writing tag_cig_out file')
    end if


    status=0
    return	

  end function RotCig_tilt


!subroutine fft_1d_data_par(fft_id,sign,n1,i_thread,data_par)
!type (fft_type) fft_id
!integer :: n1,i_thread
!integer :: sign
!complex, dimension (:,:) :: data_par
!
!call fft_1d_data(fft_id,sign,n1,data_par(:,i_thread))
!!
!end subroutine fft_1d_data_par


end module RotCig_mod
