C   The Program modifies a map
C   FILES:
C
C      5 - Terminal "LCD"
C      6 - Lineprinter "LPT"
C     10 - Input map . "MAPIN1"
c     12 - Output map.  "MAPOUT"
C
C
      PROGRAM MAPMOD
      PARAMETER (NUVLIM=160000)
      DIMENSION RHOSEC1(NUVLIM)
      CHARACTER*80 TITL1,TITLE
      COMMON /IO/ LCD,LPT,LMAPIN1,LMAPOUT
      COMMON /MAP1/NSEC1,IUVW1(3),MXYZ1(3),NW11,
     &   NU11,NU21,NV11,NV21,CELL1(6),LSPGRP1,LMODE1
     &  ,RHMIN1,RHMAX1
      CHARACTER MAPIN1*6,MAPOUT*6
      DATA MAPIN1/'MAPIN1'/,MAPOUT/'MAPOUT'/
C
C   DEFINE LOGICAL UNITS
      LCD=5
      LPT=6
      LMAPIN1=10
      LMAPOUT=12
      READ(LCD,100) TITLE
 100  FORMAT(A80)
      WRITE(LPT,110) TITLE
 110  FORMAT(/10X,'***        MAP MODIFICATION       ***',
     $ /,20X,'VERSION  MAY/93'////30X,A80//)
C
	CALL MRDHDR(LMAPIN1,MAPIN1,TITL1,NSEC1,IUVW1,MXYZ1,NW11,NU11,
     &   NU21,NV11,NV21,CELL1,LSPGRP1,LMODE1,RHMIN1,RHMAX1,RHMEAN1,
     &   rmsmp1)
C
C number of GRID POINTS in the FAST AXIS is NU
C number of GRID POINTS in the MEDIUM AXIS is NV
      NU=MXYZ1(IUVW1(1))
      NV=MXYZ1(IUVW1(2))
      NW=MXYZ1(IUVW1(3))
      call mwrhdr(LMAPOUT,title,nsec1,iuvw1,mxyz1,nw11,nu11,nu21,
     &  nv11,nv21,cell1,lspgrp1,lmode1)
C
      CALL MODMAP(RHOSEC1,NU,NV,NW)
      STOP
C
      END
C*******************************************************
      SUBROUTINE MODMAP(RHOSEC1,NU,NV,NW)
      PARAMETER (NULIM=400,NVLIM=400)
      DIMENSION RHOSEC1(NU,NV)
      COMMON /IO/ LCD,LPT,LMAPIN1,LMAPOUT
      COMMON /MAP1/NSEC1,IUVW1(3),MXYZ1(3),NW11,
     &   NU11,NU21,NV11,NV21,CELL1(6),LSPGRP1,LMODE1
     &  ,RHMIN1,RHMAX1
C...
      IF(NU.GT.NULIM .OR. NV.GT.NVLIM)GOTO 1200
      READ(LCD,*)BOT,POWER,C
      IF(ABS(BOT).LT.0.00001)BOT=0.0
      IF(ABS(POWER).LT.0.00001)POWER=2.0
      RHOMIN=100000
      RHOMAX=-100000
      RHMAXMID=-100000
      RHOMEAN=0.0
      rhosq=0.0
      DO 600 LSEC=NW11,NSEC1+NW11-1
      RHMINSEC=100000
      RHMAXSEC=-100000
      CALL MPOSN(LMAPIN1,LSEC)
      CALL MGULPR(LMAPIN1,RHOSEC1,IER)
      IF(IER .NE. 0) GOTO 1300
      IF(ABS(C).LT.0.00001.AND.LSEC.EQ.NW11) THEN
      DO 70 J=(NV21-NV11)/4,(NV21-NV11)*3/4
      DO 80 I=(NU21-NU11)/4,(NU21-NU11)*3/4
      IF(RHOSEC1(I,J) .GT. RHMAXMID)RHMAXMID=RHOSEC1(I,J)
 80   CONTINUE
 70   CONTINUE
      C=RHMAXMID*RHMAX1
      END IF
      DO 700 J=1,NV21-NV11+1
      DO 800 I=1,NU21-NU11+1
      IF(RHOSEC1(I,J).LT.BOT)RHOSEC1(I,J)=BOT
      RHOSEC1(I,J)=rhosec1(i,j)**POWER
      IF(RHOSEC1(I,J).GT.C)RHOSEC1(I,J)=C
      RHOMEAN=RHOMEAN+RHOSEC1(I,J)
      RHOSQ=RHOSQ+RHOSEC1(I,J)*RHOSEC1(I,J)
      IF(RHOSEC1(I,J) .LT. RHMINSEC)RHMINSEC=RHOSEC1(I,J)
      IF(RHOSEC1(I,J) .GT. RHMAXSEC)RHMAXSEC=RHOSEC1(I,J)
800   CONTINUE
700   CONTINUE
      CALL MSPEW(LMAPOUT,RHOSEC1)
      IF(RHMINSEC.LT.RHOMIN)RHOMIN=RHMINSEC
      IF(RHMAXSEC.GT.RHOMAX)RHOMAX=RHMAXSEC
600   CONTINUE
      WRITE(LPT,1100)BOT,POWER,C
1100  FORMAT(' LOWER LIMIT =',F10.5,',   POWER =',F6.3,
     &   ',   UPPER LIMIT =',F15.5)
      call mclose(LMAPOUT,RHOMIN,RHOMAX,rhomean,rhosq)
      RETURN
1200  STOP ' *** ARRAY SIZE NOT BIG ENOUGH ***'
1300  STOP ' *** MAP READING ERROR ***'
      END