      SUBROUTINE ll2imjm(BETA,beta_ETA,slimit,elimit,valnull)
C-----------------------------------------------------------------------
      INCLUDE "parmeta"
      INCLUDE "wave.inc"
      PARAMETER (IMT=2*IM-1,JMT=JM/2+1,LM1=LM-1,LP1=LM+1)
      PARAMETER (IMJM=IM*JM-JM/2)
      PARAMETER (DTR = .01745329)
      INCLUDE "MAPOT.comm"
C-----------------------------------------------------------------------
      DIMENSION beta     (NXP,NYP)
      DIMENSION beta_K   (IMJM)
      DIMENSION beta_ETA (IM,JM)
      DIMENSION GLAT(IM,JM), GLON(IM,JM), SM(IM,JM)
      DIMENSION GLAT1D(IMJM), GLON1D(IMJM), SM_K(IMJM)
      REAL slimit,elimit,valnull
C-----------------------------------------------------------------------
      DO J = 1, NYP
      DO I = 1, NXP
         if (beta(i,j).lt.slimit) beta(i,j) = valnull
         if (beta(i,j).gt.elimit) beta(i,j) = valnull
      ENDDO
      ENDDO
C
      DO K = 1, IMJM
         GLAT1D(K) = 0.
         GLON1D(K) = 0.
         SM_K(K) = 1.
         beta_K (K) = 0.018 
      ENDDO  
C               
      DO J = 1, JM
      DO I = 1, IM
         beta_ETA(I,J) = 0.018
      ENDDO
      ENDDO
C     
      TPH0=TPH0D*DTR
      TLM0=TLM0D*DTR
      CTPH0=COS(TPH0)
      STPH0=SIN(TPH0)
      TPHD=SBD-DPHD
      DLMD2=2.*DLMD
      DO J=1,JM
         KHL=IM*(J-1)-(J-1)/2+1
         KHH=IM* J   -    J/2
         TPHD=TPHD+DPHD
         DO K=KHL,KHH
         TLMD=WBD+DLMD2*(K-KHL)*MOD(J,2)+(DLMD+(K-KHL)*DLMD2)*MOD(J+1,2)
            TLMR=TLMD*DTR
            TPHDR=TPHD*DTR
            CALL RTLLNEW(TLMR,TPHDR,TLM0D,DTR,CTPH0,STPH0,ALMD,APHD)
            GLON1D(K)=ALMD
            GLAT1D(K)=APHD
            X=ALMD-RWEST
            Y=APHD-RSOUT
            IN=X/DELTAL+1
            JN=Y/DELTAF+1
            X=X-(IN-1)*DELTAL
            Y=Y-(JN-1)*DELTAF
            P=X/DELTAL
            Q=Y/DELTAF
            PQ=P*Q
            I00=IN
            I10=IN+1
            I01=IN
            I11=IN+1
            J00=JN
            J10=JN
            J01=JN+1
            J11=JN+1
C----------------- BILINEAR INTERPOLATION -------------------
             IF(IN.ge.1 .AND. IN.le.(NXP-1) .AND.
     &          JN.ge.1 .AND. JN.le.(NYP-1) ) THEN
c            IF(GLON1D(K).ge.RWEST .AND. GLON1D(K).le.REAST .AND.
c     &         GLAT1D(K).ge.RSOUT .AND. GLAT1D(K).le.RNORT .AND.
c     &         SM_K(K) .eq. 1.) THEN      
                beta_K(K) = beta(IN,JN)
c               beta_K(K) = beta(I00,J00)
c     *                     +P*(beta(I10,J10)-beta(I00,J00))
c     *                     +Q*(beta(I01,J01)-beta(I00,J00))
c     *                     +PQ*(beta(I00,J00)-beta(I10,J10)
c     *                     -beta(I01,J01)+beta(I11,J11))
cpk               if (k.ge.50000 .AND. k.le.50010) then
cpk                print*,'k,in,jn,beta_K,beta,GLON1D(K),GLAT1D(K)'
cpk                print*,
cpk     .   k,in,jn,beta_K(K),beta(IN,JN),GLON1D(K),GLAT1D(K)
cpk               endif
            ENDIF
         ENDDO
      ENDDO
C------------------------------------------------------------
      CALL CONH12(beta_K(1),beta_ETA(1,1))
      CALL CONH12(GLAT1D, GLAT)
      CALL CONH12(GLON1D, GLON)
C------------------------------------------------------------
cpk      bmin = 100000.
cpk      bmax = -100000.
cpk      bmin0 = 100000.
cpk      bmax0 = -100000.
cpk      DO J = 1, JM
cpk      DO I = 1, IM
cpk         if (beta_ETA(i,j) .lt. bmin) bmin = beta_ETA(i,j)
cpk         if (beta_ETA(i,j) .gt. bmax) bmax = beta_ETA(i,j)
cpk      ENDDO
cpk      ENDDO
cpk      DO J0 = 1, NYP
cpk      DO I0 = 1, NXP
cpk         if (beta(i0,j0) .lt. bmin0) bmin0 = beta(i0,j0)
cpk         if (beta(i0,j0) .gt. bmax0) bmax0 = beta(i0,j0)
cpk      ENDDO
cpk      ENDDO
cpk      print*,'bmin-bmax-bmin0-bmax0 :',bmin,bmax,bmin0,bmax0
cpk      print*,'------------------------------------------------------'
C------------------------------------------------------------
      end
C&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
      SUBROUTINE CONH12(HI1,HI2)
C     ******************************************************************
C     *                                                                *
C     *  ROUTINE FOR REORDERING THIBUE-HEIGHT-POINT-1-DIMENSIONAL      *
C     *  MATRICES FOR 2-DIMENSIONAL INDEXING                           *
C     *                                                                *
C     ******************************************************************
      INCLUDE "parmeta"
      INCLUDE "wave.inc"
      PARAMETER (IMT=2*IM-1,JMT=JM/2+1,LM1=LM-1,LP1=LM+1)
      PARAMETER (IMJM=IM*JM-JM/2)
      INCLUDE "MAPOT.comm"
C-----------------------------------------------------------------------
                             D I M E N S I O N
     & HI2    (IM,JM),HI1    (IMJM)
C-----------------------------------------------------------------------
      K=0
              DO 100 J=1,JM
          DO 100 I=1,IM-1+MOD(J,2)
      K=K+1
  100  HI2(I,J)=HI1(K)
              DO 101 J=2,JM-1,2
  101  HI2(IM,J)=HI2(IM-1,J)
C-----------------------------------------------------------------------
                             RETURN
                             END
C&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
      SUBROUTINE RTLLNEW (TLMD,TPHD,TLM0D,DTR,CTPH0,STPH0,ALMD,APHD)
      STLM=SIN(TLMD)
      CTLM=COS(TLMD)
      STPH=SIN(TPHD)
      CTPH=COS(TPHD)
!
      SPH=CTPH0*STPH+STPH0*CTPH*CTLM
      APH=ASIN(SPH)
      APHD=APH/DTR
      ANUM=CTPH*STLM
      DENOM=(CTLM*CTPH-STPH0*SPH)/CTPH0
      RELM=ATAN2(ANUM,DENOM)
      ALMD=RELM/DTR+TLM0D
!
      IF(ALMD.GT. 180.)    ALMD=ALMD-360.
      IF(ALMD.LT.-180.)    ALMD=ALMD+360.
!
      RETURN
      END
