c$DEBUG
c----------------------------------------------------------------------------
c
c   QPRINC - QUADRUPOLE TENSOR IN PRINCIPAL INERTIAL SYSTEM FROM
c            PRINCIPAL QUADRUPOLE COMPONENTS AND ANGLES BETWEEN THE TWO
c            AXIS SYSTEMS
c
C
C       ver. 23a.V.1997                           ----- Zbigniew KISIEL -----
C                          __________________________________________________
C                         | Institute of Physics, Polish Academy of Sciences |
C                         | Al.Lotnikow 32/46, Warszawa, POLAND              |
C                         |                             kisiel@ifpan.edu.pl  |
C                         |     http://info.ifpan.edu.pl/~kisiel/prospe.htm  |
C_________________________/--------------------------------------------------
C
c
      IMPLICIT REAL*8 (A-H,O-Z)
      con=0.017453293D0
c
      WRITE(*,3344)
3344  FORMAT(1X//'  ',76(1H_)/' |',T79,'|'/
     * ' |  QPRINC - QUADRUPOLE CONVERSION FROM ',
     * 'PRINCIPAL QUADRUPOLE AXES TO',T79,'|'/
     * ' |                                      ',
     * 'PRINCIPAL INERTIAL AXES',T79,'|'/
     * ' |',76(1H_),'|'/'  version 23a.V.1997',T64,'Zbigniew KISIEL'/)
c
      open(2,file='QPRINC.INP',status='OLD',err=6)
c
c  Data file:       line1:   chiz chix chiy
c                   line2:   theta.za theta.zb theta.zc
c                   line3:   theta.xa theta.xb theta.xc
c                   line4:   theta.ya theta.yb theta.yc
c
8     read(2,*,err=8,end=5)chiz,chix,chiy
      read(2,*)tza,tzb,tzc
      read(2,*)txa,txb,txc
      read(2,*)tya,tyb,tyc
c
      cza=cos(tza*con)
      czb=cos(tzb*con)
      czc=cos(tzc*con)
      cxa=cos(txa*con)
      cxb=cos(txb*con)
      cxc=cos(txc*con)
      cya=cos(tya*con)
      cyb=cos(tyb*con)
      cyc=cos(tyc*con)
c     
      chia=chix*cxa**2+chiy*cya**2+chiz*cza**2
      chib=chix*cxb**2+chiy*cyb**2+chiz*czb**2
      chic=chix*cxc**2+chiy*cyc**2+chiz*czc**2
c
      chiab=chix*cxa*cxb+chiy*cya*cyb+chiz*cza*czb
      chiac=chix*cxa*cxc+chiy*cya*cyc+chiz*cza*czc
      chibc=chix*cxb*cxc+chiy*cyb*cyc+chiz*czb*czc
c
      write(*,'(1x/1x,a)')
     *  'Diagonal components of principal quadrupole tensor:' 
      write(*,3)chiz,chix,chiy
3     format(1x/' chiz, chix, chiy =  ',3f12.5)
      write(*,4)tza,tzb,tzc,chia,chiab,chiac,
     1          txa,txb,txc,chiab,chib,chibc,
     1          tya,tyb,tyc,chiac,chibc,chic
4     format(1x/' Rotation matrix',23x,
     * ' Inertial quadrupole tensor'/
     * /' /',3f10.5,' \      /',3f12.5,' \'/
     *  ' |',3f10.5,' |  ->  |',3f12.5,' |'/
     *  ' \',3f10.5,' /      \',3f12.5,' /'/)
      write(*,10)1.5d0*chia,0.25d0*(chib-chic)
10    format(1x,t41,'      (3/2)chi.aa =',f12.5/
     *       1x,t41,'(chi.bb-chi.cc)/4 =',f12.5)
      goto 7
c
6     write(*,'(1x//'' **** ERROR opening QPRINC.INP'',a//)')char(7)
      goto 7
5     write(*,'(1x//'' **** EOF in QPRINC.INP found before any data''
     *  ,a//)')char(7)
c
7     stop
      end
c
c----------------------------------------------------------------------------
c----------------------------------------------------------------------------