text:       SUBROUTINE BSYNBplatt(NALLIN)
  *
  *-----------------------------------------------------------------------
  *
  * BSYNB is merely a continuation of BSYN, data transfered
  * via scratch-file 14
  *
  * B takes line-absorption coefficients generated by A and
  * calculates the spectra.
  *
  * Export version  1988-03-24  ********* Olof Morell *** Uppsala
  *
  *-----------------------------------------------------------------------
  *
        include 'spectrum.inc'
  *
        INTEGER    MAXNL
        PARAMETER (MAXNL=200)
        INTEGER    NLBL
        PARAMETER (NLBL=100)
  *
        CHARACTER*20  LELE(100)
        character*50 mcode
        DIMENSION ION(100),SRXLAL(2*MAXNL),
       &          GFELOG(100),ABUND(100),ETAK(2*MAXNL),ETAR(ndp*2*MAXNL),
       &          ETAD(ndp),PROF(100),XC(ndp),CHIE(100)
        real fluxme(nlbl),icenter(nlbl)
        real mum,prf,iprf
        doubleprecision XL1,XL2,eqwidth
        doubleprecision XLBEG,XLEND,DEL,XLB(100),XL(100),XLAL(2*MAXNL)
        COMMON/ATMOS/ T(NDP),PE(NDP),PG(NDP),XI(NDP),MUM(NDP),RO(NDP),
       &  nnNTAU
  *
        COMMON /CTRAN/X(NDP),S(NDP),BPLAN(NDP),XJ(NDP),HFLUX(NDP),XK(NDP)
       & ,FJ(NDP),SOURCE(NDP),TAUS(NDP),DTAUS(NDP),JTAU0,JTAU1,ISCAT
        COMMON/CSURF/ HSURF,Y1(NRAYS)
        COMMON/CANGLE/ NMY,XMY(nrays),XMY2(nrays),WMY(nrays)
        COMMON/TAUC/ TAU(ndp),DTAULN(ndp),NTAU
        COMMON/PIECES/ XL1,XL2,DEL,EPS,NMX,NLBLDU,IINT,XMYC,IWEAK
  *
  * extension for large number of wavelengths and lines (monster II)
        doubleprecision xlambda
        common/large/ xlambda(lpoint),maxlam,ABSO(NDP,lpoint),
       & absos(ndp,lpoint),absocont(ndp,lpoint)
  *
        dimension fcfc(lpoint),y1cy1c(lpoint),
       & xlm(lpoint),jlcont(lpoint)
  
        logical findtau1,hydrovelo
        real velocity
        common/velo/velocity(ndp),hydrovelo
        logical debug
        data debug /.false./ 
  *
  * NLBLDU is a dummy, real val. of NLBL is set as parameter
  *
        DATA eqwidth/0./,profold/0./
  
        PI=3.141593
  *
  * Initiate angle quadrature points and number of wavelengths
  * per block (NLBL)
  *
        NMY=NMX
        CALL GAUSI(NMY,0.,1.,WMY,XMY)
        DO 1 I=1,NMY
          XMY2(I)=XMY(I)*XMY(I)
      1 CONTINUE
  *
  * Initiate mode of calculation
  * IINT =1  : intensity at MY=XMYC
  * IINT =0  : flux
  * XL1      : wavelength where syntetic spectrum starts
  * XL2      : wavelength where syntetic spectrum stops
  * DEL      : wavelength step
  * IWEAK =1 : weak line approximation for L/KAPPA le. EPS
  * note  XL2.gt.XL1
  *
        iweak=0
        IF(IINT.GT.0) THEN
          NMY=NMY+1
          XMY(NMY)=XMYC
          wmy(nmy)=0.
          WRITE(7,200) XMYC,XL1,XL2,DEL
        END IF
        IF(IINT.EQ.0) WRITE(7,201) XL1,XL2,DEL
        IF(IWEAK.GT.0) WRITE(7,202) EPS
  *
  * Read model atmosphere
  *
        READ(14,rec=1) MCODE,nlcont,xlm(1),BPLAN,XC,S,XI
  cc       print*,'bsynbplatt read 14 ',mcode,nlcont,xlm(1),xi
        WRITE(7,203) MCODE(1:lenstr(mcode))
  *
  *
        if (iint.gt.0) then
  * header for intensity file
  *
          write(66) mcode
          write(66) 1.,1.
          write(66) maxlam
        endif
  
  *
  * Continuum calculations:
  *
        do 1963 jc=1,nlcont
          READ(14,rec=jc) MCODE,idum,xlm(jc),BPLAN,XC,S,XI
          if (debug) then
            print*,'bsynbplatt read 14 ',mcode,jc,idum,xlm(jc)
            print*,' XC     S     BPlan'
            do k=1,ntau
              print*,xc(k),s(k),bplan(k)
            enddo
          endif
        DO 9 K=1,NTAU
          X(K)=XC(K)
      9 CONTINUE
         
        call traneqplatt(0)
        Y1CY1C(jc)=Y1(NMY)
        FCFC(jc)=4.*HSURF*pi
        IF(IINT.LE.0) WRITE(7,204) fcFC(jc),xlm(jc)
        IF(IINT.GT.0) WRITE(7,205) Y1Cy1c(jc),xlm(jc)
        if (debug) then
          print*,' continuum;  lambda = ',jc,xlm(jc),fcfc(jc)
        endif
  1963  continue
  *
  * Start loop over wavelength blocks, each containing NLBL
  * wavelength points
  *
          numb=0
          DO 39 jjj=0,maxlam,nlbl
           xlend=xlambda(jjj+nlbl)
           do 390 j=1,nlbl
            xl(j)=xlambda(j+jjj)
            if (xl(j).le.xl2) then
            IF(IWEAK.LE.0.OR.IINT.LE.0) THEN
              DO 30 K=1,NTAU
  * the continuum opacity is already included in abso
  ccc              X(K)=XC(K)+ABSO(K,J+jjj)
                X(K)=ABSO(K,J+jjj)
                s(k)=absos(k,j+jjj)
                xlsingle=xl(j)
                BPLAN(k)=BPL(T(k),xlsingle)
     30       CONTINUE
  c            if (debug) then
  c              print*,' line;  lambda = ',xl(j)
  c            endif
  cc            if (abs(xl(j)-5240.41d0).lt.1.d-4.or.
  cc     &          abs(xl(j)-5242.49d0).lt.1.d-4.or.
  cc     &          abs(xl(j)-5241.62d0).lt.1.d-4)  then
  cc              print*,'lambda = ',xl(j),' calling traneqplatt',idebug
  cc              idebug=1
  cc            else
                idebug=0
  cc            endif
              call traneqplatt(idebug)
  * interpolate continuum flux
              do 3691 jc=2,nlcont
                if(xlm(jc)-xlsingle.gt.0.) then
                  jjc=jc-1
                  goto 3692
                endif
  3691        continue
  3692        continue
  	    if (nlcont.gt.1) then
                jjc=min(jjc,nlcont-1)
                fc=(fcfc(jjc+1)-fcfc(jjc))/
       &           (xlm(jjc+1)-xlm(jjc))*(xlsingle-xlm(jjc)) + fcfc(jjc)
                y1c=(y1cy1c(jjc+1)-y1cy1c(jjc))/
       &           (xlm(jjc+1)-xlm(jjc))*(xlsingle-xlm(jjc)) + y1cy1c(jjc)
                if (debug) then
                  print*,'interpolating cont: ',jjc,xlm(jjc),xlsingle,
       &                                        xlm(jjc+1)
                  print*,'fluxes ' ,fcfc(jjc),fc,fcfc(jjc+1)
                endif
  	    else
                fc=fcfc(1)
                y1c=y1cy1c(1)
              endif
              PRF=4.*pi*HSURF/FC
  * starting with version 12.1, flux is not divided by pi anymore.
  * F_lambda integrated over lambda is sigma.Teff^4
  *            fluxme(j)=hsurf*4.
              fluxme(j)=hsurf*4.*pi
              IF(IINT.GT.0) then
                iPRF=Y1(NMY)/Y1C
                icenter(j)=y1(nmy)
              endif
              PROF(J)=1.-PRF
            ELSE
              DO 31 K=1,NTAU
                ETAD(K)=ABSO(K,J+jjj)
                IF(ETAD(K).GT.EPS) GOTO 32
     31       CONTINUE
              CALL TRANW(NTAU,TAU,XMYC,BPLAN,XC,ETAD,DELI)
              PROF(J)=DELI/Y1C
              GOTO 39
     32       DO 33 K=1,NTAU
                X(K)=ABSO(K,J+jjj)
                s(k)=absos(k,j+jjj)
                xlsingle=xl(j)
                BPLAN(k)=BPL(T(k),xlsingle)
     33       CONTINUE
  c            if (debug) then
  c              print*,' line;  lambda = ',xl(j)
  c            endif
               print*,'debug ',xlsingle
              call traneqplatt(0)
  * interpolate continuum flux
              do 3693 jc=2,nlcont
                if(xlm(jc)-xlsingle.gt.0.) then
                  jjc=jc-1
                  goto 3694
                endif
  3693        continue
  3694        continue
  	    if (nlcont.gt.1) then
                jjc=min(jjc,nlcont-1)
                fc=(fcfc(jjc+1)-fcfc(jjc))/
       &           (xlm(jjc+1)-xlm(jjc))*(xlsingle-xlm(jjc)) + fcfc(jjc)
                y1c=(y1cy1c(jjc+1)-y1cy1c(jjc))/
       &           (xlm(jjc+1)-xlm(jjc))*(xlsingle-xlm(jjc)) + y1cy1c(jjc)
              else
                fc=fcfc(1)
                y1c=y1cy1c(1)
  	    endif
              PRF=4.*pi*HSURF/FC
  * starting with version 12.1, flux is not divided by pi anymore.
  * F_lambda integrated over lambda is sigma.Teff^4
  *            fluxme(j)=hsurf*4.
              fluxme(j)=hsurf*4.*pi
              IF(IINT.GT.0) then
                iPRF=Y1(NMY)/Y1C
                icenter(j)=y1(nmy)
              endif
              PROF(J)=1.-PRF
            END IF
  ccc          eqwidth=eqwidth+(prof(j)+profold)*del/2.
  ccc          profold=prof(j)
  *
  * find depth where tau_lambda=1
            if (hydrovelo) then
              findtau1=.true.
            else
              findtau1=.false.
            endif
            if (findtau1) then
            taulambda=tau(1)*(x(1)+s(1))
            do k=2,ntau
              taulambda=taulambda+(tau(k)-tau(k-1))*(x(k)+s(k)+
       &                   x(k-1)+s(k-1))*0.5
              if (taulambda.ge.1.) then
                 print333,xl(j),k,tau(k),taulambda,t(k),ro(k)
  333            format(f10.3,x,i3,x,1pe11.4,x,1pe11.4,x,0pf7.1,x,
       &                1pe11.4,0p)
  
                 goto 1966
              endif
            enddo
  1966      continue
            endif
            endif
  *
  390      continue
  cc         PRINT 2341,NUMB
  ccc        WRITE(7,206) PROF
           IF(IINT .LE. 0) DUMMY=FC
           IF(IINT .GT. 0) DUMMY=Y1C
  
  *
  * Write spectrum on file for convolution with instrument profile
  *
           do iplez=1,nlbl
             if (xl(iplez).le.xl2) then
               plezflux=1.-prof(iplez)
  
               if (iint.eq.0) then
  * fluxme is sigma teff^4
                 write(46,1964) xl(iplez),plezflux,fluxme(iplez)
  1964           format(f11.3,1x,f10.5,1x,1pe12.5)
  
  C OUTPUT IN CASE OF intensity calculation : STORE LIMB FUNCTION
               else
  * We add intensity at center of disk in spectrum output.
                 write(46,1965) xl(iplez),plezflux,fluxme(iplez),
       &                        icenter(iplez)
  1965           format(f11.3,1x,f10.5,2(1x,1pe12.5))
                 write(66) xl(iplez),nmy,fluxme(iplez),
       &                (xmy(nlx),y1(nlx),nlx=1,nmy)
               endif
  
             endif
           enddo
  *
  * Check whether there are more wavelength blocks
  *
           NUMB=NUMB+1
  
  39      CONTINUE
  *
  * End loop over wavelength blocks
  *
  ccc        WRITE(7,207) NALLIN
  *
  *
  ccc      DO 442 LIN=1,NALLIN
  ccc        WRITE(7,208) LELE(LIN),ION(LIN),XLB(LIN),CHIE(LIN),GFELOG(LIN),
  ccc     &               ABUND(LIN)
  ccc  442 CONTINUE
  ccc      print 1111,eqwidth*1000.
  1111    format(' eqwidth: ',f9.1,' mA')
          call clock
        RETURN
  *
    100 FORMAT(4X,I1,6X,I3)
    207 FORMAT(' SPECTRUM CALCULATED FOR THE FOLLOWING ',I3,' LINES'
       &      /' ELEMENT LAMBDA XI(EV) LOG(GF) LOG(ABUND) LINE NO')
    208 FORMAT('   ',A2,I2,F9.2,F5.2,1X,F6.2,3X,F7.2,4X,I5)
    200 FORMAT(' INTENSITY SPECTRUM AT MU=',F6.2,' BETWEEN',F10.3,
       &      ' AND',F10.3,' A WITH STEP',F6.3,' A')
    201 FORMAT(' FLUX SPECTRUM BETWEEN',F10.3,' A AND ',F10.3,' A WITH',
       &       ' STEP',F6.3,' A')
    202 FORMAT(' ***WEAK LINE APPROX FOR L/KAPPA L.T.',F7.4)
    203 FORMAT(' MODEL ATMOSPHERE:',A)
    204 FORMAT(' CONTINUUM FLUX=',E12.4, ' AT ',f10.2,'AA')
    205 FORMAT(' CONTINUUM INTENSITY=',E12.4,' at ',f10.2,'AA')
    206 FORMAT(1X,10F8.4)
   2341 FORMAT(1X,'Block number: ',I8)
  *
        END
  