SUBROUTINE initialize_synthetic_camera(params,F)
IMPLICIT NONE
TYPE(KORC_PARAMS), INTENT(IN) :: params
TYPE(FIELDS), INTENT(IN) :: F
! REAL(rp) :: aperture ! Aperture of the camera (diameter of lens)
! in meters
REAL(rp) :: start_at
!! in Seconds
REAL(rp) :: Riw
!! Radial position of inner wall
INTEGER, DIMENSION(2) :: num_pixels
!! Number of pixels (X,Y)
REAL(rp), DIMENSION(2) :: sensor_size
!! (horizontal,vertical)
REAL(rp) :: focal_length
REAL(rp), DIMENSION(2) :: position
!! Position of camera (R,Z)
REAL(rp) :: incline
REAL(rp) :: lambda_min
!! Minimum wavelength in cm
REAL(rp) :: lambda_max
!! Maximum wavelength in cm
INTEGER :: Nlambda
LOGICAL :: camera_on, integrated_opt, toroidal_sections, photon_count
INTEGER :: ntor_sections
REAL(rp) :: xmin, xmax, ymin, ymax, DX, DY
INTEGER :: ii
NAMELIST /SyntheticCamera/ camera_on,Riw,num_pixels,sensor_size,&
focal_length,&
position,incline,lambda_min,lambda_max,Nlambda,integrated_opt,&
toroidal_sections,ntor_sections,start_at,photon_count
open(unit=default_unit_open,file=TRIM(params%path_to_inputs),&
status='OLD',form='formatted')
read(default_unit_open,nml=SyntheticCamera)
close(default_unit_open)
! write(*,nml=SyntheticCamera)
cam%camera_on = camera_on
if (cam%camera_on) then
if (params%mpi_params%rank .EQ. 0) then
write(6,'(/,"* * * * * * * * * * * * * * * * * *")')
write(6,'("* Initializing synthetic camera *")')
end if
! cam%aperture = aperture
cam%start_at = start_at
cam%Riw = Riw
cam%num_pixels = num_pixels
cam%sensor_size = sensor_size
cam%pixel_area = PRODUCT(sensor_size/num_pixels)
cam%focal_length = focal_length
cam%position = position
cam%incline = C_PI*incline/180.0_rp
cam%horizontal_angle_view = ATAN2(0.5_rp*cam%sensor_size(1),&
cam%focal_length)
cam%vertical_angle_view = ATAN2(0.5_rp*cam%sensor_size(2),&
cam%focal_length)
cam%lambda_min = lambda_min ! In meters
cam%lambda_max = lambda_max ! In meters
cam%Nlambda = Nlambda
cam%Dlambda = (cam%lambda_max - cam%lambda_min)/REAL(cam%Nlambda,rp)
ALLOCATE(cam%lambda(cam%Nlambda))
cam%photon_count = photon_count
cam%integrated_opt = integrated_opt
cam%toroidal_sections = toroidal_sections
if (cam%toroidal_sections) then
cam%ntor_sections = ntor_sections
else
cam%ntor_sections = 1_idef
end if
cam%r = (/COS(0.5_rp*C_PI - cam%incline),-SIN(0.5_rp*C_PI - &
cam%incline),0.0_rp/)
do ii=1_idef,cam%Nlambda
cam%lambda(ii) = cam%lambda_min + REAL(ii-1_idef,rp)*cam%Dlambda
end do
ALLOCATE(cam%pixels_nodes_x(cam%num_pixels(1)))
ALLOCATE(cam%pixels_nodes_y(cam%num_pixels(2)))
ALLOCATE(cam%pixels_edges_x(cam%num_pixels(1) + 1))
ALLOCATE(cam%pixels_edges_y(cam%num_pixels(2) + 1))
xmin = -0.5_rp*cam%sensor_size(1)
xmax = 0.5_rp*cam%sensor_size(1)
DX = cam%sensor_size(1)/REAL(cam%num_pixels(1),rp)
do ii=1_idef,cam%num_pixels(1)
cam%pixels_nodes_x(ii) = xmin + 0.5_rp*DX + REAL(ii-1_idef,rp)*DX
end do
do ii=1_idef,cam%num_pixels(1)+1_idef
cam%pixels_edges_x(ii) = xmin + REAL(ii-1_idef,rp)*DX
end do
ymin = -0.5_rp*cam%sensor_size(2)
ymax = 0.5_rp*cam%sensor_size(2)
DY = cam%sensor_size(2)/REAL(cam%num_pixels(2),rp)
do ii=1_idef,cam%num_pixels(2)
cam%pixels_nodes_y(ii) = ymin + 0.5_rp*DY + REAL(ii-1_idef,rp)*DY
end do
do ii=1_idef,cam%num_pixels(2)+1_idef
cam%pixels_edges_y(ii) = ymin + REAL(ii-1_idef,rp)*DY
end do
! Initialize ang variables
ALLOCATE(ang%eta(cam%num_pixels(1)))
! angle between main line of sight and a given pixel along the x-axis
ALLOCATE(ang%beta(cam%num_pixels(1)))
ALLOCATE(ang%psi(cam%num_pixels(2)+1_idef))
! angle between main line of sight and a given pixel along the y-axis
do ii=1_idef,cam%num_pixels(1)
ang%eta(ii) = ABS(ATAN2(cam%pixels_nodes_x(ii),cam%focal_length))
if (cam%pixels_edges_x(ii) .LT. 0.0_rp) then
ang%beta(ii) = 0.5_rp*C_PI - cam%incline - ang%eta(ii)
else
ang%beta(ii) = 0.5_rp*C_PI - cam%incline + ang%eta(ii)
end if
end do
do ii=1_idef,cam%num_pixels(2)+1_idef
ang%psi(ii) = ATAN2(cam%pixels_edges_y(ii),cam%focal_length)
end do
ang%threshold_angle = ATAN2(cam%Riw,-cam%position(1))
ang%threshold_radius = SQRT(cam%Riw**2 + cam%position(1)**2)
ALLOCATE(ang%ax(cam%num_pixels(1)))
! angle between main line of sight and a given pixel along the x-axis
ALLOCATE(ang%ay(cam%num_pixels(2)))
! angle between main line of sight and a given pixel along the y-axis
do ii=1_idef,cam%num_pixels(1)
ang%ax(ii) = ATAN2(cam%pixels_nodes_x(ii),cam%focal_length)
end do
ang%min_ax = MINVAL(ang%ax)
ang%max_ax = MAXVAL(ang%ax)
do ii=1_idef,cam%num_pixels(2)
ang%ay(ii) = ATAN2(cam%pixels_nodes_y(ii),cam%focal_length)
end do
ang%min_ay = MINVAL(ang%ay)
ang%max_ay = MAXVAL(ang%ay)
if (cam%incline.GT.0.5_rp*C_PI) then
ang%ac = cam%incline - 0.5_rp*C_PI
else
ang%ac = 0.5_rp*C_PI - cam%incline
end if
! Initialize poloidal plane parameters
pplane%grid_dims = num_pixels
ALLOCATE(pplane%nodes_R(pplane%grid_dims(1)))
ALLOCATE(pplane%nodes_Z(pplane%grid_dims(2)))
! * * * * * * * ALL IN METERS * * * * * * *
IF (TRIM(params%field_model(1:10)) .EQ. 'ANALYTICAL') THEN
pplane%Rmin = F%Ro - F%AB%a
pplane%Rmax = F%Ro + F%AB%a
pplane%Zmin = -F%AB%a
pplane%Zmax = F%AB%a
ELSE
pplane%Rmin = MINVAL(F%X%R)
pplane%Rmax = MAXVAL(F%X%R)
pplane%Zmin = MINVAL(F%X%Z)
pplane%Zmax = MAXVAL(F%X%Z)
END IF
pplane%DR = (pplane%Rmax - pplane%Rmin)/REAL(pplane%grid_dims(1),rp)
pplane%DZ = (pplane%Zmax - pplane%Zmin)/REAL(pplane%grid_dims(2),rp)
do ii=1_idef,pplane%grid_dims(1)
pplane%nodes_R(ii) = pplane%Rmin + &
0.5_rp*pplane%DR + REAL(ii-1_idef,rp)*pplane%DR
end do
do ii=1_idef,pplane%grid_dims(2)
pplane%nodes_Z(ii) = pplane%Zmin + &
0.5_rp*pplane%DZ + REAL(ii-1_idef,rp)*pplane%DZ
end do
! * * * * * * * ALL IN METERS * * * * * * *
call save_synthetic_camera_params(params)
if (params%mpi_params%rank .EQ. 0) then
write(6,'("* Synthetic camera ready! *")')
write(6,'("* * * * * * * * * * * * * * * * * *",/)')
end if
end if
END SUBROUTINE initialize_synthetic_camera