SUBROUTINE sumdem(itype, switch, ear, ne, abun, dens, z, & ninputt, inputt, dem, ifl, qtherm, & velocity, photar, status) INTEGER itype, switch, ne, ninputt, ifl, status REAL ear(0:*), abun(*), inputt(*), dem(*), photar(*) REAL dens, z, velocity LOGICAL qtherm c Subroutine to calculate the summed emission from plasma with a c DEM given by the input array. c Arguments : c itype I i: type of plasma emission file c 1 = R-S c 2 = Mekal c 3 = Meka c 4 = APEC c switch I i: 0==calculate, 1==interpolate, 2=apec interpolate c ear R i: model energy ranges c ne I i: number of model energies c abun R i: abundances c dens R i: density (cm^-3) c z R i: redshift c ninputt I i: number of temperatures c inputt R i: temperatures c dem R i: emission measures for input temperatures c ifl I i: dataset number (unused at present) c qtherm R i: apply thermal broadening (APEC models only) c velocityR i: gaussian velocity broadening (APEC models only) c photar R r: spectrum c status I r: 0==OK c Entry DEMPLT CHARACTER cmd(*)*(*) REAL y(*) INTEGER iery(*) INTEGER npts, nvec, ncmd C Entry DEMSIZ INTEGER nsize c Local variables INCLUDE '../../inc/xspec.inc' INTEGER TOTEL, MKEL, RSEL, APEL PARAMETER (TOTEL=30, MKEL=14, RSEL=12, APEL=13) REAL andgrv(TOTEL), nabun(TOTEL) INTEGER mkelt(MKEL), rselt(RSEL), apelt(APEL) INTEGER i INTEGER ndem, itdem, iedem CHARACTER*2 celt(TOTEL) CHARACTER contxt*127 CHARACTER*20 name(4), pname REAL fgabnd EXTERNAL fgabnd SAVE ndem, itdem, iedem DATA celt/'H ', 'He', 'Li', 'Be', 'B ', 'C ', 'N ', 'O ', 'F ', & 'Ne', 'Na', 'Mg', 'Al', 'Si', 'P ', 'S ', 'Cl', 'Ar', & 'K ', 'Ca', 'Sc', 'Ti', 'V ', 'Cr', 'Mn', 'Fe', 'Co', & 'Ni', 'Cu', 'Zn'/ DATA andgrv/1.00e+00, 9.77e-02, 1.45e-11, 1.41e-11, 3.98e-10, & 3.63e-04, 1.12e-04, 8.51e-04, 3.63e-08, 1.23e-04, & 2.14e-06, 3.80e-05, 2.95e-06, 3.55e-05, 2.82e-07, & 1.62e-05, 3.16e-07, 3.63e-06, 1.32e-07, 2.29e-06, & 1.26e-09, 9.77e-08, 1.00e-08, 4.68e-07, 2.45e-07, & 4.68e-05, 8.32e-08, 1.78e-06, 1.62e-08, 3.98e-08/ DATA rselt /2, 6, 7, 8, 10, 12, 14, 16, 18, 20, 26, 28/ DATA apelt /2, 6, 7, 8, 10, 12, 13, 14, 16, 18, 20, 26, 28/ DATA mkelt /2, 6, 7, 8, 10, 11, 12, 13, 14, 16, 18, 20, 26, 28/ DATA ndem /0/ DATA name /'Raymond-Smith', 'Mewe-Kaastra-Liedahl', & 'Mewe-Kaastra', 'APEC'/ status = 0 pname = name(itype) c Check for invalid switch and itype combinations IF ( (switch .EQ. 0 .AND. itype .EQ. 1) .OR. & (switch .EQ. 0 .AND. itype .EQ. 4) .OR. & (switch .EQ. 1 .AND. itype .EQ. 3) ) THEN CALL xwrite('SUMDEM: This error should not occur', 2) RETURN ENDIF c These two lines are required to avoid problems running on c DUnix 4.0 for -O2 or higher ! WRITE(contxt,'(a,i2)') 'SUMDEM : itype = ', itype CALL xwrite(contxt,25) c set abundances - we use the Anders & Grevesse abundance c scale internally so shift from the Solar abundance table in use. IF ( itype .EQ. 1 ) THEN DO i = 1, RSEL nabun(i) = abun(i) * fgabnd(celt(rselt(i))) & / andgrv(rselt(i)) ENDDO ELSEIF ( itype .EQ. 4 ) THEN DO i = 1, APEL nabun(i) = abun(i) * fgabnd(celt(apelt(i))) & / andgrv(apelt(i)) ENDDO ELSE DO i = 1, MKEL nabun(i) = abun(i) * fgabnd(celt(mkelt(i))) & / andgrv(mkelt(i)) ENDDO ENDIF C If itype =2 or 3 but switch =2 then we are using one of the mekal C models but want to apply the apec code. Need to correct the nabun C array for the extra element included in mekal IF ( (itype .EQ. 2 .OR. itype .EQ. 3) .AND. switch .EQ. 2 ) THEN DO i = 6, APEL nabun(i) = nabun(i+1) ENDDO ENDIF C If the calculate option is set... IF ( switch .EQ. 0 ) THEN CALL smdem0(itype, ear, ne, nabun, dens, z, ninputt, inputt, & dem, ifl, photar, status) C Else if the interpolate option is set... ELSEIF ( switch .EQ. 1 ) THEN CALL smdem1(itype, ear, ne, nabun, dens, z, ninputt, inputt, & dem, ifl, photar, status) C Else if the APED interpolate option is set... ELSEIF ( switch .EQ. 2 ) THEN CALL smdem2(itype, ear, ne, nabun, dens, z, ninputt, inputt, & dem, ifl, qtherm, velocity, photar, status) ELSE CALL xwrite('SUMDEM : Unknown value of switch', 5) ENDIF IF ( ndem .NE. ninputt ) THEN ndem = ninputt CALL udmget(ndem, 6, itdem, status) contxt = 'Failed to get memory for tdem' IF ( status .NE. 0 ) GOTO 999 CALL udmget(ndem, 6, iedem, status) contxt = 'Failed to get memory for edem' IF ( status .NE. 0 ) GOTO 999 ENDIF DO i = 1, ninputt MEMR(itdem+i-1) = inputt(i) MEMR(iedem+i-1) = dem(i) ENDDO 999 CONTINUE IF ( status .NE. 0 ) CALL xwrite(contxt, 10) RETURN C ******** ENTRY demsiz(nsize) nsize = ndem*2 RETURN C ******** ENTRY demplt(y, iery, npts, nvec, cmd, ncmd) DO i = 1, ndem y(i) = MEMR(itdem+i-1) y(i+ndem) = MEMR(iedem+i-1) ENDDO nvec = 2 npts = ndem iery(1) = 0 iery(2) = 0 cmd(1) = 'LA X Temperature (keV)' cmd(2) = 'LA Y Emission Measure' cmd(3) = 'MARK 4 ON 2' cmd(4) = 'LINE ON 2' cmd(5) = 'LSTYLE 2 ON 2' cmd(6) = 'LOG X' cmd(7) = 'LA F '//pname ncmd = 7 RETURN END