* PACKAGE GTLAPLA   !" 球面調和関数フィルタ
*********************************************************************
      PROGRAM GTLAPLA
*
#include         "gtsinc.F"
#include         "gzsize.F"
#include         "gpainc.F"              !" /GPAWRK/
      COMMON     /GMWORK/ MWORK
      REAL       MWORK  ( IJKDIM )
*
      CHARACTER  HHEAD( NDC )*(NCC)
      REAL       GDATA( IJKDIM )
*
      CHARACTER  HFILE( 1 ) *(NFILN)
      DATA       HFILE / 'gtool.out' /
      DATA       IFILE / 50 /
      DATA       JFILE / 60 /
*
      CHARACTER  OUT    *(NFILN)
      DATA       OUT    / 'gtool.out' /
      CHARACTER  ITEM   *(NCC)
      CHARACTER  UNIT   *(NCC)
      CHARACTER  TITLE  *(NCC*2)
      CHARACTER  DSET   *(NCC)
      CHARACTER  EDIT   *(NCC)
      CHARACTER  ETTL   *(NCC)
      DATA       ITEM, UNIT, TITLE, DSET, EDIT, ETTL /6*' '/
      LOGICAL    GRESET
      DATA       GRESET / .FALSE. /
      LOGICAL    HELP
      DATA       HELP   / .FALSE. /
*
      EXTERNAL   LAPL
*     
      NAMELIST  /OPTION/ OUT,
     &                   ITEM, UNIT, TITLE, DSET, EDIT, ETTL, GRESET,
     &                   HELP, HFILE
*
      CALL OPTARG ( 91, 'OPTION', 'HFILE', NOPT, NFILE )
      READ (91,OPTION,IOSTAT=IOS)
      CLOSE(91)
      IF ( IOS.NE.0 .OR. HELP ) THEN
         WRITE(6,OPTION)
         STOP
      ENDIF
*
      CALL GURNTF ( HFILE( 1 ), OUT  , '$GTTMPDIR/gtool.in' )
*
      CALL GTOPEN
      CALL GTSIZE ( HHEAD, IJKDIM )
      CALL GMSIZE ( IJKDIM  )
*
      CALL GFROPN ( IFILE , HFILE( 1 ) )
      CALL GFWOPN ( JFILE ,  OUT )
*
      IL=LENC(OUT)
      WRITE (0,*) 'output='//OUT(1:IL)
*
 1100 CONTINUE
         CALL   GFREAD
     O        ( HHEAD , GDATA , IEOD  ,
     I          IFILE , 1               )
*
         IF ( IEOD .EQ.0 ) THEN
*
            CALL GMCAL1
     I         ( LAPL  ,
     M           HHEAD , GDATA ,
     I           EDIT  , ETTL  )
*
            IF ( ITEM .NE. ' ' ) THEN
               CALL GHCSET( HHEAD , 'ITEM', ITEM )
            ENDIF
            IF ( UNIT .NE. ' ' ) THEN
               CALL GHCSET( HHEAD , 'UNIT', UNIT )
            ENDIF
            IF ( TITLE .NE. ' ' ) THEN
               CALL GHCSTS( HHEAD , 'TITL', TITLE )
            ENDIF
            IF ( DSET .NE. ' ' ) THEN
               CALL GHCSET( HHEAD , 'DSET', DSET )
            ENDIF
            IF ( GRESET ) THEN
               CALL GHRSGP( HHEAD  )
            ENDIF
*
            CALL  GFWRIT
     I                 ( HHEAD , GDATA ,
     I                   JFILE , 1     , 0       )
*
      GOTO 1100
         ENDIF
*
      STOP
      END
*********************************************************************
      SUBROUTINE LAPL
     I         ( HH    , GDATA , 
     O           HHF   , GDATAF,
     I           IMAX  , JMAX  , KMAX  )
*
      CHARACTER  HH     ( * )*(*)        !" ヘッダー(入力)
      REAL       GDATA  ( * )            !" データ(入力)
      CHARACTER  HHF    ( * )*(*)        !" ヘッダー(出力)
      REAL       GDATAF ( * )            !" データ(出力)
*
      PARAMETER  ( IDIM=3, JDIM=2, KDIM=1 )
#include         "gpainc.F"              !" /GPAWRK/
*"    COMMON /GPAWK?/ WDATA , ZDATA , ZWORK ,
*"                    PNM   , DPNM  , TRIGS , IFAX  ,
*"                    NMO   , GW    , COSLAT,
*"                    QPNM  , QDPNM , QGW   , QSINLA
*"    REAL ER / 6370.E3 /
*
*"         < 球関数の準備 >
*
      CALL GPASPS
     I         ( IMAX  , JMAX  , KMAX  , ER    ,
     O           PNM   , DPNM  , TRIGS , IFAX  ,
     O           NMO   , GW    , COSLAT,
     O           LMAX  , MMAX  , NMAX  , MINT  ,
     O           NMDIM , JMXHF ,
     W           QPNM  , QDPNM , QGW   , QSINLA         )
*
*"         < スペクトルに >
*
      CALL SPG2W
     M         ( WDATA ,
     I           GDATA ,
     C           PNM   , NMO   , TRIGS , IFAX  , GW    ,
     F           '    ', 'POS' ,
     D           IMAX  , JMAX  , KMAX  , IMAX  , JMAX  ,
     D           LMAX  , MMAX  , NMAX  , MINT  , NMDIM , JMXHF ,
     W           ZDATA , ZWORK                                   )
*
      DO 1100 K = 1, KMAX
         DO 1100 N = 0, NMAX
            POW = 0.
            DO 1200 M = 0, MIN(MMAX,N)
               NMI = (N-M)*(MMAX+1)*2 + M*2 
               WDATA(NMO(1+NMI)+(K-1)*NMDIM) 
     &       = WDATA(NMO(1+NMI)+(K-1)*NMDIM)*N*(N+1.)
               WDATA(NMO(2+NMI)+(K-1)*NMDIM)
     &       = WDATA(NMO(2+NMI)+(K-1)*NMDIM)*N*(N+1.)
 1200       CONTINUE 
 1100 CONTINUE 
*
*"         < 格子点値に戻す >
*
      CALL SPW2G
     M         ( GDATAF,
     I           WDATA ,
     C           PNM   , NMO   , TRIGS , IFAX  ,
     F           '    ', 'POS' ,
     D           IMAX  , JMAX  , KMAX  , IMAX  , JMAX  ,
     D           LMAX  , MMAX  , NMAX  , MINT  , NMDIM , JMXHF , 
     W           ZDATA , ZWORK                                   )
*
      RETURN
      END