

! Sorting and averaging of file .abn

      dimension ABUNDI(1000,10),TAUBAR(1000,10)
      DIMENSION EPOT(1000),WDTH(1000),vts(10),ni(1000,10),gland(1000)
	character*4 id,el(100000),elo(600),id_m,str*80,fn*80
	real*8 wave
        real wbden, cbden, glden
        logical lwbden, lcbden, lglden



	  m=iargc()          ! Take file name & files opening
	   call getarg(1,fn)
	  k=index(fn,'.')
	  fn=fn(1:k)//'abn'
	   open(1,file=fn,status='old')
	   open(3,file='temp',form='unformatted',recl=80,access='direct',
     &	   status='unknown')
	  fn=fn(1:k)//'res'
	   open(2,file=fn,status='unknown')

	   nr=1
	   read(1,*)
	   read(1,*)

1	read(1,'(a80)',end=2)str   ! How many lines ?
	write(3,rec=nr)str
	read(str,'(1x,a4)')id
	el(nr)=id
	nr=nr+1
	go to 1

2	nr=nr-1     ! Start sorting
	  id_m=el(1)

10	nlines=0        ! Select lines for one element
	nvt=0
	 do i=1,nr
	  if(el(i).ne.id_m)go to 11
	  read(3,rec=i)str
	read(str,3)id,wave,e,vt,gf,glan,w,tau,a
3	format(1x,a4,f11.4,f9.4,f5.2,f8.3,4f7.2)
	  if(nvt.eq.0)then
	   nvt=1
	    ivt=1
	     vts(ivt)=vt
	   endif
	  do k=1,nvt
	   if(vt.eq.vts(k))go to 12
	  enddo
	   nvt=nvt+1
	   ivt=nvt
	   vts(ivt)=vt
	  go to 13
12	   ivt=k
13	   if(ivt.eq.1)then
   	    nlines=nlines+1	      
            EPOT(NLINES)=E
            WDTH(NLINES)=W
	    gland(nlines)=glan
	   endif
          ABUNDI(NLINES,IVT)=A
          TAUBAR(NLINES,IVT)=TAU
	  ni(nlines,ivt)=i
	  el(i)='    '
11	continue
	enddo

	write(2,'(30x,a4,/)')id_m  ! Output & statisticks
	  do ivt=1,nvt
	   write(2,"(20x,'Vturb=',f5.1,' km/sec',/)")vts(ivt)
	      SUMM=0.
	      CSUMXY=0.
	      CSUMX=0.
	      CSUMXX=0.
	      WSUMXY=0.
	      WSUMX=0.
	      WSUMXX=0.
	      glsumxx=0.
	      glsumxy=0.
	      glsumx=0.

	    do i=1,nlines
	     nrc=ni(i,ivt)
	     read(3,rec=nrc)str
	     read(str,3)id,wave,e,vt,gf,glan,w,tau,a
	      if(vt.ne.vts(ivt))go to 15
	     write(2,14)wave,e,gf,glan,w,tau,a
14	format(f11.4,f9.4,f8.3,f7.3,3f7.2)
	      if(nlines.eq.1)then
	             DEV=0.
        	    WRITE(2,33)NLINES,A,DEV
	            go to 15
	      endif
	      SUMM=SUMM+ABUNDI(I,IVT)
	      CSUMXY=CSUMXY+ABUNDI(I,IVT)*EPOT(I)
	      CSUMX=CSUMX+EPOT(I)
	      CSUMXX=CSUMXX+EPOT(I)**2
	      WSUMXY=WSUMXY+ABUNDI(I,IVT)*WDTH(I)
	      WSUMX=WSUMX+WDTH(I)
	      WSUMXX=WSUMXX+WDTH(I)**2
	      glsumxy=glsumxy+abundi(i,ivt)*gland(i)
	      glsumx=glsumx+gland(i)
	      glsumxx=glsumxx+gland(i)*gland(i)				
	    enddo
	      X=NLINES
              WBDEN = X*WSUMXX-WSUMX**2
              IF (ABS(WBDEN) .GT. 1E-9) THEN
        	WB=(X*WSUMXY-WSUMX*SUMM) / WBDEN
                LWBDEN = .TRUE.
              ELSE
                LWBDEN = .FALSE.
              END IF
              CBDEN = X*CSUMXX-CSUMX**2
              IF (ABS(CBDEN) .GT. 1E-9) THEN
                CB=(X*CSUMXY-CSUMX*SUMM) / CBDEN
                LCBDEN = .TRUE.
              ELSE
                LCBDEN = .FALSE.
              END IF
              GLDEN = x*glsumxx-glsumx**2
              IF (ABS(GLDEN) .GT. 1E-9 .AND. GLSUMXX .NE. 0.) THEN
                glb = (x*glsumxy-glsumx*summ) / GLDEN
                LGLDEN = .TRUE.
	      ELSE
                LGLDEN = .FALSE.
              END IF
	      AVABN=SUMM/X
	      SUM=0.
	      DO 32 I=1,NLINES
32 	       SUM=SUM+(AVABN-ABUNDI(I,IVT))**2
	      DEV=SQRT(SUM/X)
	      WRITE(2,33)NLINES,AVABN,DEV
   33 FORMAT(//26H ***** The abundance from I3,9H lines is,F7.2,
     1 3H+/-,F5.2)
      IF (LCBDEN) THEN
        WRITE(2,6009)CB
 6009   FORMAT(5X,' Slope=',1pE10.3,' for abundance vs. potential')
      ELSE
        WRITE (2,*) 'Cannot calculate abundance vs. potential'
      END IF
      IF (LWBDEN) THEN
        WRITE(2,6008)WB
 6008   FORMAT(5X,' Slope=',1pE10.3,' for abundance vs. widths')
      ELSE
        WRITE (2,*) 'Cannot calculate abundance vs. widths'
      END IF
      IF (LGLDEN) THEN
	write(2,6010)glb
6010	format(5x,' Slope=',1pe10.3,' for abundance vs. Lande factor',
     *  /)  
      ELSE
        WRITE (2,*) 'Cannot calculate abundance vs. Lande factor'   
      END IF

15	continue
	enddo
	
	  do i=1,nr      ! Next element
	    if(el(i).ne.'   ')then
	     id_m=el(i)
	     go to 10
	    endif
	  enddo 		   

	close(3,status='delete')	
      END

