*23456789a123456789b123456789c123456789d123456789e123456789f123456789g12
*
*
      program math
      parameter (ntemp=3)
      character spectrum(ntemp)*6
      dimension temp(ntemp)
      data temp/5000.,6000.,7000./
      data spectrum/'spec01','spec02','spec03'/
*
      h=6.6260755e-34
      c=2.99792458e+8
      bolt=1.380658e-23
*
      xcon1=h*c/bolt 
      xcon=xcon1/(1.e+4*1.e-6) 
      print*,'xcon=',xcon
*
      print*
      con=5.e+0
      x=5.e+0
  110 continue
      xold=x 
      x=con*(1.-exp(-x))
      print920,xold,x
  920 format(2f20.16)
      if(abs(x-xold) .ne. 0.) go to 110
      xmax=x
*
      print*
      wcon=h*c/(bolt*1.e-6*xmax) 
      print*,'wcon=',wcon
*
      pi=acos(-1.)
      fcon=pi*2.*h*c**2    ! this turns it into mks flux
      do 410 i=1,ntemp
      open(unit=1,file=spectrum(i),status='unknown')
      write(1,'(a4)') 'END:'
      do 420 w=0.1,1.0001,.001
      ww=w*1.e-6   ! conversion to meters
      xx=xcon1/(ww*temp(i)) 
      ff=( (fcon/ww**5)/(exp(xx)-1.) )*1.e-6 ! To make it micron integrable
      write(1,930) w,ff
  930 format(2e15.7)
  420 continue
      close(unit=1)
  410 continue
*
      end
*23456789a123456789b123456789c123456789d123456789e123456789f123456789g12
