33#if defined(HAVE_NETCDF)
54 integer,
parameter :: &
55 INTEGRATE_NONE = -1, &
68 subroutine pes_out_velocity_map(pesK, file, namespace, space, Lk, ll, how, spacing, coord_system, pmesh)
69 real(real64),
intent(in) :: pesK(:,:,:)
70 character(len=*),
intent(in) :: file
71 type(namespace_t),
intent(in) :: namespace
72 class(space_t),
intent(in) :: space
73 real(real64),
intent(in) :: Lk(:,:)
74 integer,
intent(in) :: ll(:)
75 integer(int64),
intent(in) :: how
76 real(real64),
intent(in) :: spacing(3)
77 class(coordinate_system_t),
intent(in) :: coord_system
78 real(real64),
optional,
intent(in) :: pmesh(:,:,:,:)
81 character(len=512) :: filename
83 type(cube_function_t) :: cf
88 call cube_init(cube, ll, namespace, space, spacing, coord_system)
93 if (.not.
present(pmesh))
then
96 dk(1:space%dim) = abs(lk(2, 1:space%dim) - lk(1,1:space%dim))
99#if defined(HAVE_NETCDF)
101 if (
bitand(how, option__outputformat__netcdf) /= 0)
then
102 filename = trim(file)//
".ncdf"
103 write(
message(1),
'(a)')
'Writing netcdf format file: '
112 if (
bitand(how, option__outputformat__vtk) /= 0)
then
113 filename = trim(file)//
".vtk"
114 write(
message(1),
'(a)')
'Writing vtk format file: '
117 if (
present(pmesh))
then
121 call dvtk_out_cf(filename, namespace,
'PES_vel_map', ierr, cf, cube, dk(:),&
136 type(namespace_t),
intent(in) :: namespace
137 real(real64),
intent(in) :: arpes(:,:,:)
138 character(len=*),
intent(in) :: file
139 integer,
intent(in) :: dim
140 integer,
intent(in) :: ll(:)
141 real(real64),
intent(in) :: pmesh(:,:,:,:)
142 real(real64),
intent(in) :: Ekin(:,:,:)
144 integer :: iunit, ip, ie
145 real(real64) :: dp, length, pp(1:2), pdiff(1:2)
149 iunit =
io_open(file, namespace, action=
'write')
150 write(iunit,
'(a)')
'##################################################'
151 write(iunit,
'(a1,a18,2x,a18,2x,a18,2x,a18,2x, a18,2x,a18)')
'#', &
156 write(iunit,
'(a1,a18,2x,a18,2x,a18,2x,a18,2x, a18,2x,a18)')
'#', &
162 write(iunit,
'(a)')
'##################################################'
169 pp(1:2) = pmesh(ip,1,1,1:2)
170 if (ip > 1) pdiff = pp(1:2) - pmesh(ip-1,1,1,1:2)
171 dp = norm2(pdiff(1:2))
177 write(iunit,
'(es19.12,2x,es19.12,2x,es19.12,2x,es19.12,2x,es19.12,2x,es19.12)') &
187 write(iunit,
'(es19.12,2x,es19.12,2x,es19.12,2x,es19.12,2x,es19.12,2x,es19.12)') &
211 subroutine pes_out_velocity_map_cut(namespace, pesK, file, ll, dim, pol, dir, integrate, pos, Lk, pmesh)
213 real(real64),
intent(in) :: pesk(:,:,:)
214 character(len=*),
intent(in) :: file
215 integer,
intent(in) :: ll(:)
216 integer,
intent(in) :: dim
217 real(real64),
intent(in) :: pol(3)
218 integer,
intent(in) :: dir
219 integer,
intent(in) :: integrate
220 integer,
optional,
intent(in) :: pos(3)
221 real(real64),
optional,
intent(in) :: lk(:,:)
222 real(real64),
optional,
intent(in) :: pmesh(:,:,:,:)
224 integer :: ii, ix, iy, iunit,ldir(2), icut(3)
225 real(real64) :: kk(3),temp
226 integer,
allocatable :: idx(:,:)
227 real(real64),
allocatable :: lk_(:,:)
228 real(real64) :: rotation(1:dim,1:dim)
229 logical :: aligned_axis
231 real(real64) :: k, kkk(3), theta, phi, dphi, dk(3)
234 integer :: idone, ntodo
236 real(real64),
allocatable :: cube_f(:)
241 iunit =
io_open(file, namespace, action=
'write')
245 assert(
size(pesk, 1) == ll(1))
247 if (.not.
present(pmesh))
then
250 safe_allocate(idx(1:maxval(ll(:)), 1:3))
251 safe_allocate(lk_(1:maxval(ll(:)), 1:3))
254 dk(1:dim) = abs(lk(2,1:dim)-lk(1,1:dim))
258 call sort(lk_(1:ll(ii), ii), idx(1:ll(ii), ii))
262 aligned_axis = sum((pol-(/0 ,0 ,1/))**2) <=
m_epsilon .or. &
263 sum((pol-(/0 ,1 ,0/))**2) <=
m_epsilon .or. &
266 if (
present(pos))
then
269 icut(1:3) = ll(1:3)/2 + 1
272 if (aligned_axis .and. integrate == integrate_none)
then
285 if (
present(pmesh))
then
286 do ix = 1, ll(ldir(1))
287 do iy = 1, ll(ldir(2))
291 temp = pesk(icut(dir), ix, iy)
292 kk(1) = pmesh(icut(dir), ix, iy, ldir(1))
293 kk(2) = pmesh(icut(dir), ix, iy, ldir(2))
295 temp = pesk(ix, icut(dir), iy)
296 kk(1) = pmesh(ix, icut(dir), iy, ldir(1))
297 kk(2) = pmesh(ix, icut(dir), iy, ldir(2))
299 temp = pesk(ix, iy, icut(dir))
300 kk(1) = pmesh(ix, iy, icut(dir), ldir(1))
301 kk(2) = pmesh(ix, iy, icut(dir), ldir(2))
304 write(iunit,
'(es19.12,2x,es19.12,2x,es19.12)') &
313 do ix = 1, ll(ldir(1))
314 kk(1) = lk_(ix, ldir(1))
315 do iy = 1, ll(ldir(2))
316 kk(2) = lk_(iy, ldir(2))
320 temp = pesk(idx(icut(dir), 1), idx(ix, 2), idx(iy, 3))
322 temp = pesk(idx(ix, 1), idx(icut(dir), 2), idx(iy, 3))
324 temp = pesk(idx(ix, 1), idx(iy, 2), idx(icut(dir), 3))
327 write(iunit,
'(es19.12,2x,es19.12,2x,es19.12)') &
341 print *,
"Rotate z-axis over the zenith axis"
342 print *,rotation(1,:)
343 print *,rotation(2,:)
344 print *,rotation(3,:)
349 ntodo = product(ll(1:2))
360 if (
present(pmesh))
then
361 kk(2) = pmesh(ix, 1, 1, 1)
362 kk(3) = pmesh(1, iy, 1, 2)
369 if (
present(pmesh))
then
370 kk(1) = pmesh(ix, 1, 1, 1)
371 kk(3) = pmesh(1, iy, 1, 2)
379 if (
present(pmesh))
then
380 kk(1) = pmesh(ix, 1, 1, 1)
381 kk(2) = pmesh(1, iy, 1, 2)
391 select case (integrate)
401 theta =
atan2(norm2(kk(1:2)), kk(3))
403 kkk(1) = k *
sin(theta) *
cos(phi)
404 kkk(2) = k *
sin(theta) *
sin(phi)
405 kkk(3) = k *
cos(theta)
441 write(iunit,
'(es19.12,2x,es19.12,2x,es19.12)') &
454 write(stdout,
'(1x)')
463 safe_deallocate_a(idx)
464 safe_deallocate_a(lk_)
478 real(real64),
intent(in) :: pesk(:,:,:)
479 real(real64),
intent(in) :: lk(:,:)
480 integer,
intent(in) :: ll(:)
481 integer,
intent(in) :: dim
482 real(real64),
allocatable,
intent(out) :: cube_f(:)
483 type(
qshep_t),
intent(out) :: interp
484 real(real64),
optional,
intent(in) :: pmesh(:,:,:,:)
486 integer :: np, ii, ix, iy, iz
487 real(real64) :: kk(3)
488 real(real64),
allocatable :: kx(:),ky(:),kz(:)
493 call messages_write(
"Initializing Qshep interpolator. Be patient it may take a while... ")
496 np = ll(1)*ll(2)*ll(3)
499 if (dim < 2 .or. dim > 3)
then
500 message(1) =
"This interpolator works only for 2 <= dim <= 3."
504 safe_allocate(cube_f(1:np))
506 safe_allocate(kx(1:np))
507 safe_allocate(ky(1:np))
508 safe_allocate(kz(1:np))
517 if (
present(pmesh))
then
522 cube_f(ii) = pesk(ix,iy,iz)
524 kx(ii) = pmesh(ix,iy,iz,1)
525 ky(ii) = pmesh(ix,iy,iz,2)
526 kz(ii) = pmesh(ix,iy,iz,3)
542 cube_f(ii) = pesk(ix,iy,iz)
562 safe_deallocate_a(kx)
563 safe_deallocate_a(ky)
564 safe_deallocate_a(kz)
578 real(real64),
allocatable,
intent(inout) :: cube_f(:)
579 type(
qshep_t),
intent(inout) :: interp
585 safe_deallocate_a(cube_f)
This is the common interface to a sorting routine. It performs the shell algorithm,...
double sin(double __x) __attribute__((__nothrow__
double cos(double __x) __attribute__((__nothrow__
double atan2(double __y, double __x) __attribute__((__nothrow__
subroutine, public dcube_function_alloc_rs(cube, cf, in_device, force_alloc)
Allocates locally the real space grid, if PFFT library is not used. Otherwise, it assigns the PFFT re...
subroutine, public dcube_function_free_rs(cube, cf)
Deallocates the real space grid.
subroutine, public cube_init(cube, nn, namespace, space, spacing, coord_system, fft_type, fft_library, dont_optimize, nn_out, mpi_grp, need_partition, tp_enlarge, blocksize)
subroutine, public cube_end(cube)
type(debug_t), save, public debug
real(real64), parameter, public m_two
real(real64), parameter, public m_zero
real(real64), parameter, public m_pi
some mathematical constants
real(real64), parameter, public m_epsilon
real(real64), parameter, public m_one
subroutine, public dout_cf_netcdf(filename, ierr, cf, cube, space, spacing, transpose, unit, namespace)
Writes a cube_function in netcdf format.
subroutine, public io_close(iunit, grp)
integer function, public io_open(file, namespace, action, status, form, position, die, recl, grp)
This module is intended to contain "only mathematical" functions and procedures.
subroutine, public generate_rotation_matrix(R, ff, tt)
Generates a rotation matrix R to rotate a vector f to t.
subroutine, public messages_info(no_lines, iunit, verbose_limit, stress, all_nodes, namespace)
subroutine, public messages_new_line()
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
subroutine, public pes_out_velocity_map_cut(namespace, pesK, file, ll, dim, pol, dir, integrate, pos, Lk, pmesh)
subroutine pes_out_interpolator_init(namespace, pesK, Lk, ll, dim, cube_f, interp, pmesh)
Qshep interpolation helper function initialization. Generates the linearized version of pesK (cube_f)...
integer, parameter integrate_ky
subroutine, public pes_out_velocity_map(pesK, file, namespace, space, Lk, ll, how, spacing, coord_system, pmesh)
integer, parameter integrate_r
integer, parameter integrate_phi
integer, parameter integrate_theta
subroutine pes_out_interpolator_end(cube_f, interp)
Destroy the interpolation objects.
integer, parameter integrate_kx
subroutine, public pes_out_arpes_cut(namespace, arpes, file, dim, ll, pmesh, Ekin)
integer, parameter integrate_kz
subroutine, public qshep_end(interp)
This module is intended to contain "only mathematical" functions and procedures.
character(len=80) function, public str_center(s_in, l_in)
puts space around string, so that it is centered
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
character(len=20) pure function, public units_abbrev(this)
This module defines the unit system, used for input and output.
type(unit_system_t), public units_out
type(unit_t), public unit_one
some special units required for particular quantities
subroutine, public dvtk_out_cf_structured(filename, namespace, fieldname, ierr, cf_in, cube, unit, points, ascii)
Writes a mesh_function in VTK legacy format for structured grid data (i.e. curvilinear coordinates)....
subroutine, public dvtk_out_cf(filename, namespace, fieldname, ierr, cf_in, cube, spacing, unit)
Writes a cube_function in VTK legacy format see http: