OV-GSO-DC / dustem-wrapper_idl

Blame view

src/idl/dustem_planck_function.pro 1.82 KB
Edit Raw Normal View History
427f1205   Jean-Michel Glorian   version 4.2 merged 
1
2
3
4
5
 
function dustem_planck_function,temp,wave

;+
; NAME:
;	dustem_planck_function
261aca53   Annie Hughes   improved header 
6
 
;
427f1205   Jean-Michel Glorian   version 4.2 merged 
7
8
9
 
; PURPOSE: 
;	To calculate the Planck function in units of MJy/sr 
;
261aca53   Annie Hughes   improved header 
10
11
12
 
; CATEGORY:
;    DustEMWrap, Distributed, Low-Level, User convenience, Analysis
;
427f1205   Jean-Michel Glorian   version 4.2 merged 
13
14
15
16
17
 
; CALLING SEQUENCE: 
;	bbflux = dustem_planck_function(temp,wave) 
;
; INPUT PARAMETERS: 
;	WAVE   in microns
261aca53   Annie Hughes   improved header 
18
 
;	TEMP   Scalar giving the temperature of the Planck function in degree K
427f1205   Jean-Michel Glorian   version 4.2 merged 
19
 
;
261aca53   Annie Hughes   improved header 
20
21
 
; OUTPUT:
;	BBFLUX - Scalar or vector giving the Planck function at the specified
427f1205   Jean-Michel Glorian   version 4.2 merged 
22
 
;		wavelength points.
261aca53   Annie Hughes   improved header 
23
24
25
26
27
 
; COMMON BLOCKS:
;    None
;
; SIDE EFFECTS:
;    None
427f1205   Jean-Michel Glorian   version 4.2 merged 
28
29
30
31
32
 
;
; EXAMPLES:
;
; RESTRICTIONS:
;	Values less than approximately 1E-24 are truncated to 0.
261aca53   Annie Hughes   improved header 
33
 
;;
427f1205   Jean-Michel Glorian   version 4.2 merged 
34
 
; MODIFICATION HISTORY:
261aca53   Annie Hughes   improved header 
35
36
37
 
;    Written by JPB Apr-2011
;    Evolution details on the DustEMWrap gitlab.
;    See http://dustemwrap.irap.omp.eu/ for FAQ and help.
427f1205   Jean-Michel Glorian   version 4.2 merged 
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
 
;-

 On_error,2

 if ( N_elements(wave) LT 1 ) then begin
     print,'Syntax - bbflux = dustem_planck_function(temp,wave)'
     print,'wave in microns, Temp in K, bbflux in MJy/sr'
     return,0
  endif    

;  if ( N_elements( temp ) NE 1 ) then $
;      read,'Enter a blackbody temperture',temp

  Nwav=n_elements(wave)
  Ntemp=n_elements(temp)
  Nel=Nwav*Ntemp
  bbflux = fltarr(Nel)
;  bbflux = wave*0.

; Gives the blackbody flux (i.e. PI*Intensity) ergs/cm2/s/a
  w = wave / 1.e4   ;microns to cm
;  w = wave / 1.E8                              ; Angstroms to cm    
  c1 =  3.74185E-5             ;constants appropriate to cgs units.
  C2 =  1.43883
  val =  c2/w/temp               
  good = where( val LT 88, Ngood )    ;Avoid floating underflow
  if ( Ngood GT 0 ) then  $
      bbflux[ good ] =  C1 / ( w[good]^5 * ( exp( val[good])-1. ) )

; fact is 1.e4/1.e7/2.99792e10*1e20/!pi
  fact=1.06177e+06
  bbflux=bbflux*w*w*fact ;MJy/sr

;  return, bbflux*1.E-8              ; Convert to ergs/cm2/s/A
  RETURN,bbflux

  end