New Upstream Release - neat
Ready changes
Summary
Merged new upstream version: 2.4 (was: 2.3.2).
Resulting package
Built on 2023-02-26T10:33 (took 2m27s)
The resulting binary packages can be installed (if you have the apt repository enabled) by running one of:
apt install -t fresh-releases neat-dbgsymapt install -t fresh-releases neat
Lintian Result
Diff
diff --git a/Atomic-data/README b/Atomic-data/README
index dfd173a..9914567 100644
--- a/Atomic-data/README
+++ b/Atomic-data/README
@@ -7,7 +7,7 @@ SIII collisional data:
Transition probabilities from Mendoza and Zeippen (1982MNRAS.199.1025M) and collision strengths from Mendoza (1983IAUS..103..143M).
All other collisional data:
-CHIANTI 7.0 (Landi et al., 2012; ApJS, 744, 99)
+CHIANTI 9.0 (Dere et al., 1997; A&AS, 125, 149, Dere et al., 2019; ApJS, 241, 22)
Recombination data:
diff --git a/source/Ilines_levs b/Atomic-data/zones.dat
similarity index 100%
rename from source/Ilines_levs
rename to Atomic-data/zones.dat
diff --git a/Makefile b/Makefile
index 1f8c898..2165784 100644
--- a/Makefile
+++ b/Makefile
@@ -96,7 +96,6 @@ install: neat
test -e $(DESTDIR)$(PREFIX)/bin || mkdir -p $(DESTDIR)$(PREFIX)/bin
test -e $(MANDIR) || mkdir -p $(MANDIR)
install -m 644 Atomic-data/*.* $(DESTDIR)$(PREFIX)/share/neat
- install -m 644 source/Ilines_levs $(DESTDIR)$(PREFIX)/share/neat
install -m 644 utilities/complete_line_list $(DESTDIR)$(PREFIX)/share/neat
install -m 644 utilities/plot.sh $(DESTDIR)$(PREFIX)/share/neat
install -m 644 config/default.cfg $(DESTDIR)$(PREFIX)/share/neat
diff --git a/debian/changelog b/debian/changelog
index a65020a..f6ce57e 100644
--- a/debian/changelog
+++ b/debian/changelog
@@ -1,3 +1,9 @@
+neat (2.4-1) UNRELEASED; urgency=low
+
+ * New upstream release.
+
+ -- Debian Janitor <janitor@jelmer.uk> Sun, 26 Feb 2023 10:31:40 -0000
+
neat (2.3.2-2) unstable; urgency=medium
* changed autopkgtest to write text output instead of fits (Closes: #982795)
diff --git a/debian/patches/Fix-gcc10-FTBFS.patch b/debian/patches/Fix-gcc10-FTBFS.patch
index c40f944..e773433 100644
--- a/debian/patches/Fix-gcc10-FTBFS.patch
+++ b/debian/patches/Fix-gcc10-FTBFS.patch
@@ -7,10 +7,10 @@ Last-Update: 2021-02-02
Makefile | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
-diff --git a/Makefile b/Makefile
-index 1f8c898..539b73a 100644
---- a/Makefile
-+++ b/Makefile
+Index: neat.git/Makefile
+===================================================================
+--- neat.git.orig/Makefile
++++ neat.git/Makefile
@@ -49,7 +49,7 @@ ifeq ($(MESSAGES),yes)
endif
diff --git a/man/neat.1 b/man/neat.1
index 597b96d..3e49ae8 100644
--- a/man/neat.1
+++ b/man/neat.1
@@ -33,6 +33,10 @@ Valid values are:
Default: How
+.TP
+.B \-hb, \-\-hbeta\-flux [VALUE]
+Flux of H beta. If your line list does not include H beta, use this option to provide a value, so that NEAT can analyse the line list and compute abundances. If this command line option is used when H beta is already in the line list, the command line option will override the value in the line list.
+
.TP
.B \-he, \-\-helium\-data [VALUE]
The atomic data to use for He I abundances
@@ -75,6 +79,10 @@ Number of iterations. Default: 1
.B \-nelow, \-nemed, \-nehigh, \-telow, \-temed, \-tehigh [VALUE]
By default, electron densities and temperatures are calculated from the available diagnostics, but their values can be overridden with these commands. The required units are cm\-3 for densities and K for temperatures.
+.TP
+.B \-\-no\-omp
+This option disables parallelisation so that if uncertainties are being calculated, only one thread is used.
+
.TP
.B \-rp
When calculating Monte Carlo uncertainties, \fBneat\fR's default behaviour is to assume that all uncertainties are Gaussian. If \-rp is specified, it will compensate for the upward bias affecting weak lines described by Rola and Pelat (1994), assuming log normal probability distributions for weaker lines.
diff --git a/source/abundances.f90 b/source/abundances.f90
index 9fdb8cf..1673087 100644
--- a/source/abundances.f90
+++ b/source/abundances.f90
@@ -57,7 +57,7 @@ use mod_globals
!atomic data
- integer :: iion !# of ions in Ilines
+ integer :: iion !# of ions in zones.dat
integer :: maxlevs,maxtemps
type(atomic_data),dimension(24) :: atomicdata
real(kind=dp), dimension(21,14,44) :: heidata
@@ -185,7 +185,7 @@ use mod_globals
do srloop=1,subtract_recombination
if (calculate_extinction) then
- call calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction, dble(10000.),dble(1000.),weights%ha,weights%hg,weights%hd,hidata)
+ call calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction, dble(10000.),dble(1000.),hidata)
if (meanextinction .lt. 0.0 .or. isnan(meanextinction)) then
meanextinction = 0.d0
@@ -415,7 +415,7 @@ use mod_globals
!update extinction. DS 22/10/11
meanextinction=0
- call calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction, medtemp, lowdens,weights%ha,weights%hg,weights%hd,hidata)
+ call calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction, medtemp, lowdens, hidata)
if (meanextinction .lt. 0.0 .or. isnan(meanextinction)) then
meanextinction = 0.d0
@@ -1459,7 +1459,7 @@ if (ciiCELabund + ciiiCELabund + civCELabund .gt. 0.d0) CELicfC = CabundCEL/(cii
if (niiCELabund + niiiCELabund + nivCELabund + nvCELabund .gt. 0.d0) CELicfN = NabundCEL/(niiCELabund + niiiCELabund + nivCELabund + nvCELabund)
if (oiiCELabund + oiiiCELabund + oiiiIRCELabund + oivCELabund .gt. 0.d0) CELicfO = OabundCEL/(oiiCELabund + oiiiCELabund + oiiiIRCELabund + oivCELabund)
if (neiiiCELabund + neivCELabund + nevCELabund .gt. 0.d0) CELicfNe = NeabundCEL/(neiiiCELabund + neivCELabund + nevCELabund)
-if (siiCELabund + siiiCELabund + siiiIRCELabund + sivIRCELabund .gt. 0.d0) CELicfS = SabundCEL/(siiCELabund + siiiCELabund + siiiIRCELabund + sivIRCELabund)
+if (siiCELabund + siiiCELabund + siiiIRCELabund + sivIRCELabund .gt. 0.d0) CELicfS = SabundCEL/(siiCELabund + siiiCELabund + sivIRCELabund)
if (cliiCELabund + cliiiCELabund + clivCELabund .gt. 0.d0) CELicfCl = ClabundCEL/(cliiCELabund + cliiiCELabund + clivCELabund)
if (ariiiCELabund + arivCELabund + arvCELabund .gt. 0.d0) CELicfAr = ArabundCEL/(ariiiCELabund + arivCELabund + arvCELabund)
!recombination lines
@@ -1488,7 +1488,7 @@ if (CELicfNe .lt. 1.) then
endif
if (CELicfS .lt. 1.) then
CELicfS = 1.0
- SabundCEL = (siiCELabund + siiiCELabund + siiiIRCELabund + sivIRCELabund)
+ SabundCEL = (siiCELabund + siiiCELabund + sivIRCELabund)
endif
if (CELicfCl .lt. 1.) then
CELicfCl = 1.0
@@ -1924,7 +1924,7 @@ subroutine get_average_abundance(startion,endion,abundance)
implicit none
character(len=11) :: startion,endion
real(kind=dp) :: abundance
- real(kind=dp), dimension(9) :: weights, abundances ! Ilines_levs contains at most 9 lines per ion. this might need updating in the future
+ real(kind=dp), dimension(9) :: weights, abundances ! zones.dat contains at most 9 lines per ion. this might need updating in the future
integer :: i,j
!debugging
diff --git a/source/bashcompletion b/source/bashcompletion
index 0533145..7039b2a 100755
--- a/source/bashcompletion
+++ b/source/bashcompletion
@@ -9,7 +9,7 @@ _neat()
cur="${COMP_WORDS[COMP_CWORD]}"
prev="${COMP_WORDS[COMP_CWORD-1]}"
- opts="-i --input -u --uncertainties -n --n-iterations -e --extinction-law -c -nelow -nemed -nehigh -telow -temed -tehigh -he --helium-data -icf --ionisation-correction-scheme -v --verbosity -id --identify -idc --identify-confirm -rp -sr --subtract-recombination -cf --configuration-file --citation -o --output-dir -of --output-format"
+ opts="-i --input -u --uncertainties -n --n-iterations -e --extinction-law -c -nelow -nemed -nehigh -telow -temed -tehigh -he --helium-data -icf --ionisation-correction-scheme -v --verbosity -id --identify -idc --identify-confirm -rp -sr --subtract-recombination -cf --configuration-file --citation -o --output-dir -of --output-format -hb --hbeta-flux --no-omp"
COMPREPLY=( $(compgen -W "${opts}" -- ${cur}) )
@@ -48,7 +48,7 @@ _neat()
_filedir -d
return 0
;;
- -c|-nehigh|-nelow|-nemed|-n|--n-iterations|-tehigh|-telow|-temed|--citation)
+ -c|-nehigh|-nelow|-nemed|-n|--n-iterations|-tehigh|-telow|-temed|--citation|-hb|--hbeta-flux|--no-omp)
COMPREPLY=()
return 0
;;
diff --git a/source/commandline.f90 b/source/commandline.f90
index a15a211..89b8590 100644
--- a/source/commandline.f90
+++ b/source/commandline.f90
@@ -6,7 +6,7 @@ use mod_functions
contains
-subroutine readcommandline(runs,switch_ext,switch_he,switch_icf,meanextinction,diagnostics,verbosity,R,identifylines,identifyconfirm,nbins,norp,calculate_extinction,subtract_recombination,configfile,nperbin)
+subroutine readcommandline(runs,switch_ext,switch_he,switch_icf,meanextinction,diagnostics,verbosity,R,identifylines,identifyconfirm,nbins,norp,calculate_extinction,subtract_recombination,configfile,nperbin,hbetaflux,noomp)
implicit none
@@ -18,7 +18,8 @@ subroutine readcommandline(runs,switch_ext,switch_he,switch_icf,meanextinction,d
character :: switch_icf !switch for which ICF scheme to use
logical :: file_exists,identifylines,identifyconfirm,norp,calculate_extinction
type(diagnostic_array) :: diagnostics
- real(kind=dp) :: meanextinction, R
+ real(kind=dp) :: meanextinction, R, hbetaflux
+ logical :: noomp
#ifdef CO
print *,"subroutine: readcommandline"
@@ -165,6 +166,12 @@ subroutine readcommandline(runs,switch_ext,switch_he,switch_icf,meanextinction,d
call exit(104)
endif
endif
+ if ((trim(options(i))=="-hb" .or. trim(options(i))=="--hbeta-flux") .and. (i+1) .le. Narg) then
+ read (options(i+1),*) hbetaflux
+ endif
+ if (trim(options(i))=="--no-omp") then
+ noomp=.true.
+ endif
! to be fully implemented:
! -R : R (default 3.1) - only used with CCM at the moment
! -b : batch mode
diff --git a/source/extinction.f90 b/source/extinction.f90
index 36d0c78..8e94f0e 100644
--- a/source/extinction.f90
+++ b/source/extinction.f90
@@ -9,27 +9,27 @@ implicit none
contains
-subroutine calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction, temp, dens, weightha, weighthg, weighthd,hidata)
+subroutine calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction, temp, dens, hidata)
IMPLICIT NONE
type(line), dimension(:), intent(in) :: linelist
integer, dimension(3:40) :: H_Balmer
- real(kind=dp) :: c1, c2, c3, weightha, weighthg, weighthd, meanextinction
+ real(kind=dp) :: c1, c2, c3, meanextinction
real(kind=dp) :: temp, dens
real(kind=dp), dimension(:,:,:,:), allocatable :: hidata
- real(kind=dp), dimension(3:40) :: balmerratios,extinctioncoefficients,r1,r2,r3,r4,r5,r6
- integer :: temperaturestart,densitystart,temperatureend,densityend,i
+ real(kind=dp), dimension(3:25) :: balmerratios,extinctioncoefficients,r1,r2,r3,r4,r5,r6
+ integer :: temperaturestart,densitystart,temperatureend,densityend
real(kind=dp) :: temperatureinterpolation,densityinterpolation
!debugging
#ifdef CO
- print *,"subroutine: calc_extinction_coeffs. Te=",temp,", ne=",dens,", weights=",weightha, weighthg, weighthd
+ print *,"subroutine: calc_extinction_coeffs. Te=",temp,", ne=",dens,", weights=",linelist(H_Balmer(3:6))%weight
#endif
!interpolate HI emissivities to temp and dens
!first temperature
do temperatureend=1,ntemps
- if (temp .lt. temperatures(temperatureend)) then
+ if (temp .lt. hidatatemperatures(temperatureend)) then
exit
endif
enddo
@@ -43,7 +43,7 @@ subroutine calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction,
temperatureinterpolation=0.
else
temperaturestart=temperatureend-1
- temperatureinterpolation=(temp-temperatures(temperaturestart))/(temperatures(temperatureend)-temperatures(temperaturestart))
+ temperatureinterpolation=(temp-hidatatemperatures(temperaturestart))/(hidatatemperatures(temperatureend)-hidatatemperatures(temperaturestart))
endif
!then density which runs from 2 to 8 in the data file
@@ -67,10 +67,10 @@ subroutine calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction,
!now interpolate
- r1=hidata(temperaturestart,densitystart,3:40,2) / hidata(temperaturestart,densitystart,4,2)
- r2=hidata(temperaturestart,densityend,3:40,2) / hidata(temperaturestart,densityend,4,2)
- r3=hidata(temperatureend,densitystart,3:40,2) / hidata(temperatureend,densitystart,4,2)
- r4=hidata(temperatureend,densityend,3:40,2) / hidata(temperatureend,densityend,4,2)
+ r1=hidata(temperaturestart,densitystart,3:25,2) / hidata(temperaturestart,densitystart,4,2)
+ r2=hidata(temperaturestart,densityend,3:25,2) / hidata(temperaturestart,densityend,4,2)
+ r3=hidata(temperatureend,densitystart,3:25,2) / hidata(temperatureend,densitystart,4,2)
+ r4=hidata(temperatureend,densityend,3:25,2) / hidata(temperatureend,densityend,4,2)
r5=densityinterpolation*(r2-r1)+r1
r6=densityinterpolation*(r4-r3)+r3
@@ -80,8 +80,8 @@ subroutine calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction,
! calculate extinction coefficients
extinctioncoefficients=0
- where (H_Balmer(:) .gt. 0)
- extinctioncoefficients=log10( ( linelist(H_Balmer(:))%intensity / linelist(H_Balmer(4))%intensity )/ balmerratios(:) )/(-linelist(H_Balmer(:))%flambda)
+ where (H_Balmer(3:25) .gt. 0)
+ extinctioncoefficients=log10( ( linelist(H_Balmer(3:25))%intensity / linelist(H_Balmer(4))%intensity )/ balmerratios(:) )/(-linelist(H_Balmer(3:25))%flambda)
endwhere
! set negative values to zero
@@ -95,8 +95,6 @@ subroutine calc_extinction_coeffs(linelist,H_Balmer, c1, c2, c3, meanextinction,
! calculate weighted mean, Ha to Hd for now
-! meanextinction=sum(extinctioncoefficients*linelist(H_Balmer(:))%weight)/sum(linelist(H_Balmer(:))%weight)
-
if (sum(linelist(H_Balmer(3:6))%weight).gt.0) then
meanextinction=sum(extinctioncoefficients(3:6)*linelist(H_Balmer(3:6))%weight)/(sum(linelist(H_Balmer(3:6))%weight)-linelist(H_Balmer(4))%weight)
else
diff --git a/source/filereading.f90 b/source/filereading.f90
index e234308..066c070 100644
--- a/source/filereading.f90
+++ b/source/filereading.f90
@@ -359,7 +359,7 @@ subroutine read_celdata(ILs, ionlist)
Iint = 1
301 format(A11, 1X, A8, 1X, F7.2, 1X, A20,1X,A4,1X,A15)
- open(201, file=trim(PREFIX)//"/share/neat/Ilines_levs", status='old')
+ open(201, file=trim(PREFIX)//"/share/neat/zones.dat", status='old')
read (201,*) numberoflines !number of lines to read in is given at the top of the file
! allocate (ILs(numberoflines)) todo: restore
Ils%name=' '
diff --git a/source/helium.f90 b/source/helium.f90
index e5f9882..3f1df62 100644
--- a/source/helium.f90
+++ b/source/helium.f90
@@ -6,7 +6,7 @@ use mod_types
implicit none
real(kind=dp), dimension(:,:,:,:), allocatable :: heiidata
- real(kind=dp), dimension(:), allocatable :: temperatures
+ real(kind=dp), dimension(:), allocatable :: heiidatatemperatures
integer :: ntemps, ndens, nlevs
contains
@@ -601,7 +601,7 @@ nlevs=25 ! maximum number of levels shown in intrat data file
!allocate the emissivities array, dimensions are temperature, density, level 1, level 2
allocate(heiidata(ntemps, ndens, nlevs, nlevs))
-allocate(temperatures(ntemps))
+allocate(heiidatatemperatures(ntemps))
heiidata(:,:,:,:)=0.d0
@@ -611,7 +611,7 @@ do i=1,ntemps
do j=1,ndens
read (357,*) invar(1:6) !line at top of each block has density, charge, temperature, case, maximum level calculated and maximum level displayed
if (j.eq.1) then
- read(invar(3),"(E9.2)") temperatures(i) ! get the temperatures for later use in interpolating
+ read(invar(3),"(E9.2)") heiidatatemperatures(i) ! get the temperatures for later use in interpolating
endif
do k=nlevs,1,-1
do l=1,k-1
@@ -662,7 +662,7 @@ allocate(emissivityarray(size(heiidata,3),size(heiidata,4)))
!get the box it's in in te vs. ne grid:
do i=1,ntemps-1
- if (medtemp .ge. temperatures(i) .and. medtemp .le. temperatures(i+1)) then
+ if (medtemp .ge. heiidatatemperatures(i) .and. medtemp .le. heiidatatemperatures(i+1)) then
exit
endif
enddo
@@ -675,9 +675,9 @@ d=floor(log10(density))
if (i .lt. 1) then
t1=1
t2=1
- elseif (i .eq. size(temperatures)) then
- t1=size(temperatures)
- t2=size(temperatures)
+ elseif (i .eq. size(heiidatatemperatures)) then
+ t1=size(heiidatatemperatures)
+ t2=size(heiidatatemperatures)
else
t1=i
t2=i+1
@@ -694,8 +694,8 @@ d=floor(log10(density))
d2=d
endif
- x1=temperatures(t1)
- x2=temperatures(t2)
+ x1=heiidatatemperatures(t1)
+ x2=heiidatatemperatures(t2)
y1=real(d1+1)
y2=real(d2+1)
y=log10(density)
diff --git a/source/hydrogen.f90 b/source/hydrogen.f90
index 1dcb468..4069514 100644
--- a/source/hydrogen.f90
+++ b/source/hydrogen.f90
@@ -6,7 +6,7 @@ use mod_globals
implicit none
real(kind=dp), dimension(:,:,:,:), allocatable :: hidata
-real(kind=dp), dimension(:), allocatable :: temperatures
+real(kind=dp), dimension(:), allocatable :: hidatatemperatures
integer :: ntemps, ndens, nlevs
contains
@@ -31,7 +31,7 @@ nlevs=25 ! maximum number of levels shown in intrat data file
!allocate the emissivities array, dimensions are temperature, density, level 1, level 2
allocate(hidata(ntemps, ndens, nlevs, nlevs))
-allocate(temperatures(ntemps))
+allocate(hidatatemperatures(ntemps))
hidata(:,:,:,:)=0.d0
@@ -41,7 +41,7 @@ do i=1,ntemps
do j=1,ndens
read (357,*) invar(1:6) !line at top of each block has density, charge, temperature, case, maximum level calculated and maximum level displayed
if (j.eq.1) then
- read(invar(3),"(E9.2)") temperatures(i) ! get the temperatures for later use in interpolating
+ read(invar(3),"(E9.2)") hidatatemperatures(i) ! get the temperatures for later use in interpolating
endif
do k=nlevs,1,-1
do l=1,k-1
@@ -94,7 +94,7 @@ allocate(searcharray(size(hidata,2),size(hidata,3),size(hidata,4)))
!now find the temperature
do i=1,ntemps
- if (medtemp .lt. temperatures(i)) then
+ if (medtemp .lt. hidatatemperatures(i)) then
exit
endif
enddo
@@ -105,10 +105,10 @@ if (i>ntemps) i=ntemps ! because in do i=1,10, final value of i is 11
!if temperature is outside data limits, it is just set to the limit
!lower limit is 500K, upper is 30000K
-if (medtemp .lt. temperatures(1) .or. medtemp .gt. temperatures(size(temperatures))) then
+if (medtemp .lt. hidatatemperatures(1) .or. medtemp .gt. hidatatemperatures(size(hidatatemperatures))) then
searcharray(:,:,:)=hidata(i,:,:,:)
else
- interp_t=(medtemp-temperatures(i-1))/(temperatures(i)-temperatures(i-1))
+ interp_t=(medtemp-hidatatemperatures(i-1))/(hidatatemperatures(i)-hidatatemperatures(i-1))
searcharray(:,:,:)=hidata(i-1,:,:,:)+interp_t*(hidata(i,:,:,:)-hidata(i,:,:,:))
endif
@@ -203,7 +203,7 @@ allocate(searcharray(size(hidata,2),size(hidata,3),size(hidata,4)))
!now find the temperature
do i=1,ntemps
- if (medtemp .lt. temperatures(i)) then
+ if (medtemp .lt. hidatatemperatures(i)) then
exit
endif
enddo
@@ -214,10 +214,10 @@ if (i>ntemps) i=ntemps ! because in do i=1,10, final value of i is 11
!if temperature is outside data limits, it is just set to the limit
!lower limit is 500K, upper is 30000K
-if (medtemp .lt. temperatures(1) .or. medtemp .gt. temperatures(size(temperatures))) then
+if (medtemp .lt. hidatatemperatures(1) .or. medtemp .gt. hidatatemperatures(size(hidatatemperatures))) then
searcharray(:,:,:)=hidata(i,:,:,:)
else
- interp_t=(medtemp-temperatures(i-1))/(temperatures(i)-temperatures(i-1))
+ interp_t=(medtemp-hidatatemperatures(i-1))/(hidatatemperatures(i)-hidatatemperatures(i-1))
searcharray(:,:,:)=hidata(i-1,:,:,:)+interp_t*(hidata(i,:,:,:)-hidata(i,:,:,:))
endif
diff --git a/source/neat.f90 b/source/neat.f90
index bd2d995..2a68fa6 100644
--- a/source/neat.f90
+++ b/source/neat.f90
@@ -62,7 +62,7 @@ program neat
!atomic data
character(len=10), dimension(24) :: ionlist !list of ion names
- integer :: iion !# of ions in Ilines
+ integer :: iion !# of ions in zones.dat
integer :: maxlevs,maxtemps
type(atomic_data),allocatable :: atomicdata(:)
real(kind=dp), dimension(:,:,:), allocatable :: heidata
@@ -99,11 +99,16 @@ program neat
!OpenMP
integer :: omp_get_num_threads
+ logical :: noomp
!diagnostic array
type(diagnostic_array) :: diagnostics
+!manually set hbeta flux
+
+ real(kind=dp) :: hbetaflux
+
! default values
runs=1
@@ -131,6 +136,8 @@ program neat
configfile=""
defaultconfigfile=trim(PREFIX)//"/share/neat/default.cfg"
outputformat="fits"
+ hbetaflux=0.d0
+ noomp=.false.
iion=24 ! number of ions for which we can calculate abundances
! todo: calculate this automatically
@@ -141,15 +148,12 @@ program neat
print *
print *,gettime(),"starting code"
- call readcommandline(runs,switch_ext,switch_he,switch_icf,meanextinction,diagnostics,verbosity,R,identifylines,identifyconfirm,nbins,norp,calculate_extinction,subtract_recombination,configfile,nperbin)
+ call readcommandline(runs,switch_ext,switch_he,switch_icf,meanextinction,diagnostics,verbosity,R,identifylines,identifyconfirm,nbins,norp,calculate_extinction,subtract_recombination,configfile,nperbin,hbetaflux,noomp)
! set up filenames. if input file is FITS, output will be written to it. Output will be put in same directory as input by default.
! todo: allow non-FITS to be requested
write (outputfilename,"(A)") filename(index(filename,"/",back=.true.)+1:len(trim(filename)))
- if (trim(outputfilename(index(outputfilename,".",.true.)+1:len(filename))).ne."fits") then
- outputfilename=trim(outputfilename)//".fits"
- endif
if (trim(outputdirectory).eq."") then
outputdirectory=trim(filename(1:index(filename,"/",.true.)))
@@ -170,6 +174,52 @@ program neat
call read_text_linelist(linelist,listlength,ncols,runs)
endif
+! add hbeta if manually specified
+ if (hbetaflux.gt.0d0) then
+ if (minval(abs(linelist(:)%wavelength - 4861.33)) .lt. 0.001) then
+! overwrite, warn
+ linelist(minloc(abs(linelist(:)%wavelength - 4861.33),1))%intensity=hbetaflux
+ print *,gettime(),"warning: manually specified Hbeta flux overrides value already in line list"
+ else
+! create array one element longer
+ listlength=listlength+1
+ allocate (linelist_original(listlength+1))
+
+! copy linelist to it
+ linelist_original = linelist
+
+! deallocate linelist, reallocate one element longer
+ deallocate(linelist)
+ allocate(linelist(listlength))
+
+! copy original lines back and deallocate
+ linelist(1:listlength)=linelist_original
+ deallocate(linelist_original)
+
+! add hbeta
+ linelist(listlength)%intensity = hbetaflux
+ linelist(listlength)%abundance = 0.D0
+ linelist(listlength)%weight = 0.d0
+ linelist(listlength)%freq=0d0
+ linelist(listlength)%wavelength=4861.33
+ linelist(listlength)%wavelength_observed=4861.33
+ linelist(listlength)%int_dered=0d0
+ linelist(listlength)%int_err=0d0
+ linelist(listlength)%blend_intensity=0.d0
+ linelist(listlength)%blend_int_err=0d0
+ linelist(listlength)%zone=' '
+ linelist(listlength)%transition=' '
+ linelist(listlength)%location=0
+ linelist(listlength)%ion='H 4 '
+ linelist(listlength)%multiplet='H4 '
+ linelist(listlength)%lowerterm='2p+ 2P* '
+ linelist(listlength)%upperterm='4d+ 2D '
+ linelist(listlength)%g1=8
+ linelist(listlength)%g2=32
+
+ endif
+ endif
+
allocate (linelist_original(listlength))
! run line identifier if required
@@ -184,6 +234,7 @@ program neat
normalise = 100.d0/linelist(minloc(abs(linelist(:)%wavelength - 4861.33),1))%intensity
else
print *,gettime(),"error: no H beta detected. No further analysis possible"
+ print *,gettime()," -hb / --hbeta-flux command line option can be used to give a value"
call exit(201)
endif
@@ -310,6 +361,9 @@ program neat
maxlevs = maxval(atomicdata%nlevs)
maxtemps = maxval(atomicdata%ntemps)
+!omp setup
+ if (noomp) call omp_set_num_threads(1)
+
!now check number of iterations. If 1, line list is fine as is. If more than one, randomize the fluxes
allocate(all_results(runs))
@@ -375,6 +429,23 @@ program neat
call write_output(runs,listlength,ncols,all_linelists,all_results,verbosity,nbins,subtract_recombination)
endif
+! deallocate
+ deallocate(all_linelists)
+ deallocate(all_results)
+ deallocate(linelist)
+ deallocate(linelist_original)
+ deallocate(hidata)
+ deallocate(heidata)
+ deallocate(heiidata)
+ deallocate(atomicdata)
+ deallocate(ciirls)
+ deallocate(niirls)
+ deallocate(oiirls)
+ deallocate(neiirls)
+ deallocate(xiiirls)
+ deallocate(hidatatemperatures)
+ deallocate(heiidatatemperatures)
+
print *
print *,gettime(),"all done"
print *
diff --git a/source/output.f90 b/source/output.f90
index 935247f..68cd6de 100644
--- a/source/output.f90
+++ b/source/output.f90
@@ -31,8 +31,8 @@ subroutine write_output(runs,listlength,ncols,all_linelists,all_results,verbosit
allocate(quantity_result(runs))
quantity_result=0d0
- open (650,file=trim(filename)//"_linelist", status='replace', access='sequential', action='write')
- open (651,file=trim(filename)//"_linelist.tex", status='replace', access='sequential', action='write')
+ open (650,file=trim(outputfilename)//"_linelist", status='replace', access='sequential', action='write')
+ open (651,file=trim(outputfilename)//"_linelist.tex", status='replace', access='sequential', action='write')
if (ncols .ge. 4) then
write (650,*) " Lambda (Rest) Ion F(line) I(line) Abundance"
@@ -157,8 +157,8 @@ subroutine write_output(runs,listlength,ncols,all_linelists,all_results,verbosit
close(651)
print *,gettime(),"Linelist written to files:"
- print *," ",trim(filename),"_linelist"
- print *," ",trim(filename),"_linelist.tex"
+ print *," ",trim(outputfilename),"_linelist"
+ print *," ",trim(outputfilename),"_linelist.tex"
!now write out the summary files and all the binned data
@@ -176,8 +176,8 @@ subroutine write_output(runs,listlength,ncols,all_linelists,all_results,verbosit
!open the files and write the headers
- open (650,file=trim(filename)//"_results", status='replace', access='sequential', action='write')
- open (651,file=trim(filename)//"_results.tex", status='replace', access='sequential', action='write')
+ open (650,file=trim(outputfilename)//"_results", status='replace', access='sequential', action='write')
+ open (651,file=trim(outputfilename)//"_results.tex", status='replace', access='sequential', action='write')
write (650,*) "NEAT (nebular empirical analysis tool)"
write (650,*) "======================================"
@@ -273,8 +273,8 @@ subroutine write_output(runs,listlength,ncols,all_linelists,all_results,verbosit
close(651)
print *,gettime(),"Results written to files:"
- print *," ",trim(filename),"_results"
- print *," ",trim(filename),"_results.tex"
+ print *," ",trim(outputfilename),"_results"
+ print *," ",trim(outputfilename),"_results.tex"
end subroutine write_output
@@ -307,6 +307,10 @@ subroutine write_fits(runs,listlength,ncols,all_linelists,all_results,nbins)
print *,"subroutine: write_fits"
#endif
+ if (trim(outputfilename(index(outputfilename,".",.true.)+1:len(filename))).ne."fits") then
+ outputfilename=trim(outputfilename)//".fits"
+ endif
+
status=0
readwrite=1
varidat=0
@@ -334,18 +338,18 @@ subroutine write_fits(runs,listlength,ncols,all_linelists,all_results,nbins)
call ftibin(unit,listlength,tfields,ttype_lines,tform_lines,tunit_lines,extname,varidat,status)
- call ftpcld(unit,1,1,1,listlength,all_linelists(:,1)%wavelength,status)
- call ftpcld(unit,2,1,1,listlength,all_linelists(:,1)%wavelength_observed,status)
- call ftpcld(unit,3,1,1,listlength,all_linelists(:,1)%intensity,status)
- call ftpcld(unit,4,1,1,listlength,all_linelists(:,1)%int_err,status)
- call ftpcnd(unit,5,1,1,listlength,0.d0,status)
- call ftpcnd(unit,6,1,1,listlength,0.d0,status)
- call ftpcls(unit,7,1,1,listlength,all_linelists(:,1)%ion,status)
- call ftpcls(unit,8,1,1,listlength,all_linelists(:,1)%multiplet,status)
- call ftpcls(unit,9,1,1,listlength,all_linelists(:,1)%lowerterm,status)
- call ftpcls(unit,10,1,1,listlength,all_linelists(:,1)%upperterm,status)
- call ftpclj(unit,11,1,1,listlength,all_linelists(:,1)%g1,status)
- call ftpclj(unit,12,1,1,listlength,all_linelists(:,1)%g2,status)
+ call ftpcld(unit,1,1,1,listlength,(/all_linelists(:,1)%wavelength/),status)
+ call ftpcld(unit,2,1,1,listlength,(/all_linelists(:,1)%wavelength_observed/),status)
+ call ftpcld(unit,3,1,1,listlength,(/all_linelists(:,1)%intensity/),status)
+ call ftpcld(unit,4,1,1,listlength,(/all_linelists(:,1)%int_err/),status)
+ call ftpcnd(unit,5,1,1,listlength,(/0.d0/),status)
+ call ftpcnd(unit,6,1,1,listlength,(/0.d0/),status)
+ call ftpcls(unit,7,1,1,listlength,(/all_linelists(:,1)%ion/),status)
+ call ftpcls(unit,8,1,1,listlength,(/all_linelists(:,1)%multiplet/),status)
+ call ftpcls(unit,9,1,1,listlength,(/all_linelists(:,1)%lowerterm/),status)
+ call ftpcls(unit,10,1,1,listlength,(/all_linelists(:,1)%upperterm/),status)
+ call ftpclj(unit,11,1,1,listlength,(/all_linelists(:,1)%g1/),status)
+ call ftpclj(unit,12,1,1,listlength,(/all_linelists(:,1)%g2/),status)
else
@@ -389,7 +393,7 @@ subroutine write_fits(runs,listlength,ncols,all_linelists,all_results,nbins)
quantity_result = all_linelists(i,:)%abundance
call get_uncertainties(quantity_result, binned_quantity_result, uncertainty_array, unusual,nbins)
- call ftpcld(unit,16,i,1,1,uncertainty_array(2),status)
+ call ftpcld(unit,16,i,1,1,(/uncertainty_array(2)/),status)
call ftpcld(unit,17,i,1,1,(/uncertainty_array(2)+uncertainty_array(1)/),status)
call ftpcld(unit,18,i,1,1,(/uncertainty_array(2)-uncertainty_array(3)/),status)
enddo
@@ -433,8 +437,8 @@ subroutine write_fits(runs,listlength,ncols,all_linelists,all_results,nbins)
do j=1,166
quantity_result=resultprocessingarray(j,:)
call get_uncertainties(quantity_result, binned_quantity_result, uncertainty_array, unusual,nbins)
- call ftpcls(unit,1,j,1,1,resultprocessingtext(j,1),status)
- call ftpcld(unit,2,j,1,1,uncertainty_array(2),status)
+ call ftpcls(unit,1,j,1,1,(/resultprocessingtext(j,1)/),status)
+ call ftpcld(unit,2,j,1,1,(/uncertainty_array(2)/),status)
call ftpcld(unit,3,j,1,1,(/uncertainty_array(2)+uncertainty_array(1)/),status)
call ftpcld(unit,4,j,1,1,(/uncertainty_array(2)-uncertainty_array(3)/),status)
enddo
diff --git a/source/weights.f90 b/source/weights.f90
index 835c96b..9c59b83 100644
--- a/source/weights.f90
+++ b/source/weights.f90
@@ -26,14 +26,10 @@ subroutine setweights(configfile,weights,linelist,ILs,H_Balmer,H_Paschen,HeII_li
where (H_Balmer.gt.0)
linelist(H_Balmer)%weight = -1
- elsewhere
- linelist(H_Balmer)%weight = 0.d0
endwhere
where (H_Paschen.gt.0)
linelist(H_Paschen)%weight = -1
- elsewhere
- linelist(H_Paschen)%weight = 0.d0
endwhere
!open config file
Debdiff
[The following lists of changes regard files as different if they have different names, permissions or owners.]
Files in second set of .debs but not in first
-rw-r--r-- root/root /usr/lib/debug/.build-id/f8/1c7efbc8da70d040691da31fedc2af333a5600.debug -rw-r--r-- root/root /usr/share/neat/zones.dat
Files in first set of .debs but not in second
-rw-r--r-- root/root /usr/lib/debug/.build-id/51/b6882c1fdffdd885c9a4e9734dc718fbea1ba1.debug -rw-r--r-- root/root /usr/share/neat/Ilines_levs
No differences were encountered between the control files of package neat
Control files of package neat-dbgsym: lines which differ (wdiff format)
Build-Ids: 51b6882c1fdffdd885c9a4e9734dc718fbea1ba1 f81c7efbc8da70d040691da31fedc2af333a5600