      subroutine gray(graydat,iu,name,page1,page2,hfin,nfsm,milsec)
c                ----
c               dimensions max cartes = 1600,4000
c------------------------------------------------------------- 
      integer*4 bl(30),wh(30),zer(30),page1,page2,lb(36),
     1          ctrast(30)
      character comm*80,name*22,end*3,graydat*6,hfin*8
      dimension black(30),white(30),tab(1600,4000)
      integer*2 lec(1600)
      integer*4 head(512)
      maxdim=4000
c
1     format(a)
2     format(1x,a)
3     format(10i6)
4     format(14x,10i6)

      write(3,1)' gray'
c
      rewind(8)
      call SNAME(8,graydat,ier)
      read(8,1)comm
      read (8,3)inveri,inverj,npy,imax
        type 1,comm
        type 3,inveri,inverj,npy,imax
        if(page1.eq.0)page1=2
        if(page2.eq.0)page2=1
        if(npy.eq.0)npy=1
        npage=page1*page2*npy
      read(8,1)comm
      read(8,3)bl,wh

c                    les 30 positions recoivent des valeurs correspondant
c                    successivement a tous les tableaux rencontres dans
c                    le fichier q (ou p)
c                    ordre: Y,npol,lambda (j,i)

      do l=1,30
      ctrast(l)=0
      zer(l)=0
      enddo 


      lm=0
          do l=1,30
          black(l)=bl(l)
          white(l)=wh(l)
            if(bl(l).ne.wh(l))then
            lsup=l
            endif
          enddo
c      read(8,1)comm
c      read(8,3)ii1,ii2,jj1,jj2
c      read(8,1)comm
c      read(8,3)iz1,iz2,iz3,iz4,milsig,iprej,nt

       ii1=0
       ii2=0
       jj1=0
       jj2=0
      iz1=0
      iz2=0
      iz3=0
      iz4=0
      milsig=0
      iprej=0
      nt=0


      nplot=0
      ihv=1
      np=0
c
      nf=1
      nfm=1
          type *,'dans gray name :',name
          write(7,*)'dans gray name :',name

        rewind(iu)
c                               !!!!!!!!!!! lecture header (iu) a faire
c      im=
c      jm=
       nm=2  !!!!!!!!!!!!
c
c        type 4,im,jm,nm,kj,kh,kr
c        type *,(lb(l),l=1,lsup)

      read(iu)head!!!!!!!!!!!!!
      do i=1,30!!!!!!!!!!!!!!!
c      type *,'head(i) :',head(i)!!!!!!!!!!!!!!!!
      enddo!!!!!!!!!!!!!!!!!
      im=head(2)!!!!!!!!!!!!!!!!
      jm=head(3)!!!!!!!!!!!!!!!!!
      nm=head(4)    !!!!!!!!!!!!!!!  
      do l=1,lsup !!!!!!!!!!!!!!
      lb(l)=head(12+l)!!!!!!!!!!!!!!!
      enddo!!!!!!!!!!!!!!!!!!
10    format(' im,jm,nm,kj,kh,kr,lb(1),lb(2)... ')
      write(3,10)
        write(3,3)im,jm,nm,kj,kh,kr,(lb(l),l=1,lsup)
          if(ii1.eq.0)then
          i1=1
          i2=im
          j1=1
          j2=jm
          else
          i1=ii1
          i2=ii2
          j1=jj1
          j2=jj2
          endif
c
        llsup=lsup
        if(nm.lt.llsup)llsup=nm
c------
        do l=1,llsup
        llb=lb(l)
c
        if(bl(l).ne.wh(l))np=np+1
        if(np.gt.npy)np=1
          ihv=1

          call SUBGRAY(iu,lec,tab,l,i1,i2,im,imax,j1,j2,jm,llb,
     1          inveri,inverj,page1,page2,np,npy,milsec,
     2          black,white,ctrast,zer,
     3          iz1,iz2,iz3,iz4,milsig,iprej,nt,name,ihv,hfin,nfsm)
c

        nplot=nplot+1
        enddo

      ijtot=page1*page2*npy
      jpiv=nplot/ijtot
      jsaut=ijtot-(nplot-jpiv*ijtot)
         write(3,*)' nplot,ijtot,jsaut ',nplot,ijtot,jsaut

      if(jsaut.ne.0.and.jsaut.ne.ijtot)then
        do jst=1,jsaut
        call pgadvance
        enddo
      endif


c      if(las.ne.0)then
c         if(nplot.ne.npage.and.nf.ne.nfm)goto200
c         call pgend
c         if(las.eq.1)then
c             type *,'Print new page?'
c             pause
c             call system('lpr -Pqms gray.ps')
c             pause
c         endif        
c         if(las.eq.-1)call system('ghostview gray.ps')
c         call system('rm gray.ps')
c         if(displ.eq.4)call pgbegin(0,'gray.ps/ps',page1,page2)
c         if(displ.eq.5)call pgbegin(0,'gray.ps/vps',page1,page2)
c         np=0
c         nplot=0
c      endif

      return
      end
c-------------------------------------
      subroutine SUBGRAY(iu,lec,tab,l,i1,i2,im,imax,j1,j2,jm,llb,
     1            inveri,inverj,page1,page2,np,npy,milsec,
     2            bl,wh,ctrast,zer,
     3            iz1,iz2,iz3,iz4,milsig,iprej,nt,nom,ihv,hfin,nfsm)
      integer*2 lec(1600)
      integer*4 zer(10),nrej(90),page1,page2,ctrast(10)
      dimension tab(1600,4000),tr(6),cr(10),bl(30),wh(30)
c modif 1 mai 00  titre a 48
      character*48 titre
      character*22 nom
      character*8 hfin
      character*2 blancs
      maxdim=4000
      blancs='  '
c
      if(bl(l).eq.wh(l))goto10
      write(titre(1:22),'(a)')nom
      write(titre(23:26),'(i4)')l
      write(titre(27:32),'(i6)')llb
c        kpiv1=bl(l)
c        kpiv2=wh(l)
      write(titre(33:34),'(a)')blancs
      write(titre(35:42),'(a)')hfin
      write(titre(43:48),'(i6)')nfsm
        tr(1)=0.
        tr(2)=1.
        tr(3)=0.
        tr(4)=0.
        tr(5)=0.
        tr(6)=1.
c      do i=1,38 !!!!!!!!!!!!!!!!
c      if(titre(i:i).eq.' ')titre(i:i)='0'!!!!!!!!!!!!!
c      enddo!!!!!!!!!!!!!!
c
      ratio=(0.880/0.653)*float(page1)/float(page2)
      if(ihv.eq.1)ratio=(0.653/.880)*float(page1)/float(page2)
      is=i2-i1+1
      js=j2-j1+1
c
      ipiv=(is-1)/maxdim+1
      jpiv=(js-1)/maxdim+1
      ijpas=max0(ipiv,jpiv)
      is=(is-1)/ijpas+1
      js=(js-1)/ijpas+1
      xi1=1.
      xi2=is
      yj1=1.
      yj2=js
               piv=float(milsec)/1000.
      xxi1=0.
      xxi2=piv*(xi2-1.)*ijpas
      yyj1=0.
      yyj2=piv*(yj2-1.)*ijpas
c
      xl=0.75     
      if(npy.gt.1)xl=0.85
      if(imax.ne.0)xl=xl*float(is-1)/float(imax-1)
      yl=xl*(j2-j1)/(ratio*float(i2-i1))
c        type *, xl,yl
c
        if(npy.gt.1.and.yl.gt.0.85/float(npy))then
        piv=xl/yl
        yl=0.85/float(npy)
        xl=yl*piv
        endif
c
        if(npy.le.1.and.yl.gt.0.8)then
        piv=xl/yl
        yl=0.75
        xl=yl*piv
        endif
c
      x1=.5-0.5*xl
      x2=x1+xl
        if(npy.gt.1)then
        y1=0.04+0.90*(np-0.5)/(float(npy)+0.05)-0.5*yl
        else
        y1=0.5-0.5*yl
        endif
      y2=y1+yl
c

10    do20 j=1,jm
      read(iu)(lec(i),i=1,im)
c          --


      if(bl(l).eq.wh(l))goto20
        jtest=j-100*(j/100)
12      format(' im,jm,tab=',10i6)
c        if(jtest.eq.0)write(3,12)im,jm,(lec(i),i=1,im,200)
      if(j.lt.j1.or.j.gt.j2)goto20
      if(zer(l).ne.0)call ZERO3(iz1,iz2,iz3,iz4,milsig,iprej,nt,
     1                      lec,1,im,1,0,
     2                      ar,br,cr,item,nrjm,nrej,yy1,yy2,ifin)
       pasij=ijpas
       if(inverj.eq.0)then
        yp=(j-j1)/pasij+1
        jp=yp+0.001
        if(yp-jp.ne.0.)goto20
       else
        yp=(j2-j)/pasij+1
        jp=yp+0.001
       if(yp-jp.ne.0.)goto20
       endif
c
        do15 i=i1,i2,ijpas
        if(inveri.eq.0)then
         xp=(i-i1)/pasij+1
         ip=xp
         if(xp-ip.ne.0.)goto15
        else
         xp=(i2-i)/pasij+1
         ip=xp
         if(xp-ip.ne.0.)goto15
        endif
          tab(ip,jp)=lec(i)
          if(lec(i).eq.32000)tab(ip,jp)=32000.
15      continue
c
20    continue
c
      if(bl(l).eq.wh(l))return
      if(np.eq.1)call PGADVANCE
      call PGIDEN
      call PGVPORT(x1,x2,y1,y2)
      call PGWINDOW(xi1,xi2,yj1,yj2)
c        type *,x1,x2,y1,y2
c        type *,xi1,xi2,yj1,yj2
         piv=page2
         piv=sqrt(piv)
      call pgsch(piv)
      if(npy.ne.1)call pgsch(piv)
      if(np.eq.1.and.npy.eq.1)then
       call pgwindow(xxi1,xxi2,yyj1,yyj2)
       call PGBOX('BCInts',50.,5,'BCInts',50.,5)
       call pgwindow(xi1,xi2,yj1,yj2)
      endif
      if(np.eq.npy)call PGLABEL(' arcsec ',' ',titre)
c
      tabmax=-10000.
      tabmin=10000.
        iimin=1
        iimax=1
        jjmin=1
        jjmax=1
      do j=1,js
      do30 i=1,is
      piv=tab(i,j)
        if(abs(piv).gt.31999.)goto30
      tabmin=amin1(piv,tabmin)
      tabmax=amax1(piv,tabmax)
      if(piv.eq.tabmin)then
          iimin=i+i1-1
          jjmin=j+j1-1
      endif     
      if(piv.eq.tabmax)then
          iimax=i+i1-1
          jjmax=j+j1-1  
      endif
30    continue
31    format(' is,js,tab='2i6,1x8f6.0)
c      jtest=j-100*(j/100)
c      if(jtest.eq.0)write(3,31)is,js,(tab(i,j),i=1,is,200)
      enddo
32    format(' l,      bl,wh, imin,jmin,tabmin, imax,jmax,tabmax ='/
     1       i4,2f7.0,2i4,f7.0,2i4,f7.0)
        type 32, l,bl(l),wh(l),iimin,jjmin,tabmin,iimax,jjmax,tabmax
        write(3,32)l,bl(l),wh(l),iimin,jjmin,tabmin,iimax,jjmax,tabmax
        write(7,32) l,bl(l),wh(l),iimin,jjmin,tabmin,iimax,jjmax,tabmax

c      endif
c
      bbl=bl(l)
      wwh=wh(l)
      if(bl(l).eq.1.)then
      bbl=tabmax
      wwh=tabmin
      endif
      if(wh(l).eq.1.)then
      bbl=tabmin
      wwh=tabmax
      endif
c
      alp=(ctrast(l)+100.)/200.
      den=wwh-bbl
c      type *,alp,den
c
      do j=1,js
      do i=1,is
        if(tab(i,j).eq.32000)then
        tab(i,j)=0.5*(bbl+wwh)
        goto40
        endif
      if(alp.ne.0.5)then
      bet=(tab(i,j)-bbl)/den
        if(bet.le.alp)then
        tab(i,j)=bbl+0.5*den*bet/alp
        else
        tab(i,j)=wwh-0.5*den*(1.-bet)/(1.-alp)
        endif
      endif
40    continue
      enddo
      enddo
c
      call PGGRAY(tab,1600,4000,1,is,1,js,bbl,wwh,tr)
      return
      end
