ROHSA Program

Uses

  • program~~rohsa~~UsesGraph program~rohsa ROHSA module~mod_constants mod_constants program~rohsa->module~mod_constants module~mod_inout mod_inout program~rohsa->module~mod_inout module~mod_start mod_start program~rohsa->module~mod_start module~mod_rohsa mod_rohsa program~rohsa->module~mod_rohsa iso_fortran_env iso_fortran_env module~mod_constants->iso_fortran_env module~mod_inout->module~mod_constants module~mod_rohsa->module~mod_constants module~mod_rohsa->module~mod_start module~mod_array mod_array module~mod_rohsa->module~mod_array module~mod_functions mod_functions module~mod_rohsa->module~mod_functions module~mod_optimize mod_optimize module~mod_rohsa->module~mod_optimize module~mod_array->module~mod_constants module~mod_functions->module~mod_constants module~mod_functions->module~mod_array module~mod_functions->module~mod_optimize module~mod_optimize->module~mod_constants module~mod_optimize->module~mod_array

Calls

program~~rohsa~~CallsGraph program~rohsa ROHSA proc~read_cube read_cube program~rohsa->proc~read_cube proc~read_parameters read_parameters program~rohsa->proc~read_parameters proc~main_rohsa main_rohsa program~rohsa->proc~main_rohsa proc~read_map read_map program~rohsa->proc~read_map proc~mean_array mean_array proc~main_rohsa->proc~mean_array proc~dim2nside dim2nside proc~main_rohsa->proc~dim2nside proc~go_up_level go_up_level proc~main_rohsa->proc~go_up_level proc~set_stdmap set_stdmap proc~main_rohsa->proc~set_stdmap proc~header header proc~main_rohsa->proc~header proc~dim_data2dim_cube dim_data2dim_cube proc~main_rohsa->proc~dim_data2dim_cube proc~ender ender proc~main_rohsa->proc~ender proc~reshape_up reshape_up proc~main_rohsa->proc~reshape_up proc~mean_map mean_map proc~main_rohsa->proc~mean_map proc~std std proc~set_stdmap->proc~std proc~timestamp timestamp proc~header->proc~timestamp proc~ender->proc~timestamp

Contents

Source Code


Variables

Type AttributesNameInitial
logical :: noise

if false --> STD map computed by ROHSA with lstd and ustd (if true given by the user)

logical :: regul

if true --> activate regulation

logical :: descent

if true --> activate hierarchical descent to initiate the optimization

integer :: n_gauss

number of gaussian to fit

integer :: n_gauss_add

number of gaussian to add at each step

integer :: m

number of corrections used in the limited memory matrix by LBFGS-B

integer :: lstd

lower bound to compute the standard deviation map of the cube (if noise .eq. false)

integer :: ustd

upper bound to compute the standrad deviation map of the cube (if noise .eq. false)

integer :: iprint

print option

integer :: iprint_init

print option init

integer :: maxiter

max iteration for L-BFGS-B alogorithm

integer :: maxiter_init

max iteration for L-BFGS-B alogorithm (init mean spectrum)

real(kind=xp) :: lambda_amp

lambda for amplitude parameter

real(kind=xp) :: lambda_mu

lamnda for mean position parameter

real(kind=xp) :: lambda_sig

lambda for dispersion parameter

real(kind=xp) :: lambda_var_amp

lambda for variance amplitude parameter

real(kind=xp) :: lambda_var_mu

lambda for variance mean position parameter

real(kind=xp) :: lambda_var_sig

lambda for variance dispersion parameter

real(kind=xp) :: amp_fact_init

times max amplitude of additional Gaussian

real(kind=xp) :: sig_init

dispersion of additional Gaussian

character(len=512) :: filename_parameters

name of the parameters file (default parameters.txt)

character(len=512) :: filename

name of the data file

character(len=512) :: fileout

name of the output result

character(len=512) :: filename_noise

name of the file with STD map (if noise .eq. true)

character(len=8) :: init_option

Init ROHSA with the mean or the std spectrum

real(kind=xp), dimension(:,:,:), allocatable:: data

initial fits data

real(kind=xp), dimension(:,:), allocatable:: std_cube

standard deviation map fo the cube is given by the user


Source Code

program ROHSA

  use mod_constants
  use mod_start
  use mod_inout
  use mod_rohsa
  
  implicit none

  logical :: noise           !! if false --> STD map computed by ROHSA with lstd and ustd (if true given by the user)
  logical :: regul           !! if true --> activate regulation
  logical :: descent         !! if true --> activate hierarchical descent to initiate the optimization
  integer :: n_gauss         !! number of gaussian to fit
  integer :: n_gauss_add     !! number of gaussian to add at each step
  integer :: m               !! number of corrections used in the limited memory matrix by LBFGS-B
  integer :: lstd            !! lower bound to compute the standard deviation map of the cube (if noise .eq. false)
  integer :: ustd            !! upper bound to compute the standrad deviation map of the cube (if noise .eq. false)
  integer :: iprint          !! print option 
  integer :: iprint_init     !! print option init
  integer :: maxiter         !! max iteration for L-BFGS-B alogorithm
  integer :: maxiter_init    !! max iteration for L-BFGS-B alogorithm (init mean spectrum)
  real(xp) :: lambda_amp     !! lambda for amplitude parameter
  real(xp) :: lambda_mu      !! lamnda for mean position parameter
  real(xp) :: lambda_sig     !! lambda for dispersion parameter
  real(xp) :: lambda_var_amp !! lambda for variance amplitude parameter
  real(xp) :: lambda_var_mu  !! lambda for variance mean position parameter
  real(xp) :: lambda_var_sig !! lambda for variance dispersion parameter
  real(xp) :: amp_fact_init  !! times max amplitude of additional Gaussian
  real(xp) :: sig_init       !! dispersion of additional Gaussian

  character(len=512) :: filename_parameters !! name of the parameters file (default parameters.txt)
  character(len=512) :: filename            !! name of the data file
  character(len=512) :: fileout             !! name of the output result
  character(len=512) :: filename_noise      !! name of the file with STD map (if noise .eq. true)
  character(len=8)   :: init_option !!Init ROHSA with the mean or the std spectrum    

  real(xp), dimension(:,:,:), allocatable :: data        !! initial fits data
  real(xp), dimension(:,:), allocatable   :: std_cube    !! standard deviation map fo the cube is given by the user 

  !Print header and get filename in argument
  call get_command_argument(1, filename_parameters)
    
  !Default user parameters
  n_gauss = 6
  n_gauss_add = 0
  lambda_amp = 1._xp
  lambda_mu = 1._xp
  lambda_sig = 1._xp
  lambda_var_amp = 0._xp
  lambda_var_mu = 0._xp
  lambda_var_sig = 1._xp
  amp_fact_init = 2._xp/3._xp
  sig_init = 5._xp
  maxiter_init = 15000
  maxiter = 800
  m = 10
  noise = .false.
  regul = .true.
  descent = .false.
  lstd = 1; ustd = 20
  init_option = "mean"
  iprint = -1
  iprint_init = -1
 
  !Read parameters
  call read_parameters(filename_parameters, filename, fileout, filename_noise, n_gauss, n_gauss_add, &
       lambda_amp, lambda_mu, lambda_sig, lambda_var_amp, lambda_var_mu, lambda_var_sig, amp_fact_init, &
       sig_init, init_option, maxiter_init, maxiter, m, noise, regul, descent, lstd, ustd, iprint, iprint_init)

  !Load data
!   print*, "filename = '",trim(filename),"'"
  call read_cube(filename, data)
  
  if (noise .eqv. .true.) then
     if (filename_noise == " ") then
        print*, "--> noise = .true. (no input rms map)"
     end if
     call read_map(filename_noise, std_cube)
  end if
  
  write(*,*) ""
  write(*,*) "opening file and reading data"
  
  !Call ROHSA subroutine
  call main_rohsa(data, std_cube, fileout, n_gauss, n_gauss_add, lambda_amp, &
       lambda_mu, lambda_sig, lambda_var_amp, lambda_var_mu, lambda_var_sig, amp_fact_init, sig_init, &
       maxiter_init, maxiter, m, noise, regul, descent, lstd, ustd, init_option, iprint, iprint_init)  
   
end program ROHSA