subroutine irwpf(jatom,ion,indmol,t,u)
c     ======================================
c
c     partition functions adter Irwin (1981), ApJS. 45, 621.
c     updated with the data of Barklem & Collet (2016) 
C     set to the Irwin format by Y. Ossorio
c
c     Input: jatom - atomic number; if =0 - molecules
c            ion - ionization degree
c            indmol - index of a molecule in the new Tsuji-type
c                     indexing (from file tsuji.molec_bc2)
c            t - temperature
c     Output: u - partition function
c
c     array IRWIND(I) - the Irwin index corresponding to Tsuji
c           index I
c           if =0 - molecule I has no data in the Irwin table
c
      INCLUDE 'PARAMS.FOR'
      real*8 a(6,3,92),aa(6),am(6,500),spec(500)
      dimension irwind(478)
      save iread,a,am
c
      data irwind/
     *    0,   1,  28,   4,   2,   7,   6,   5,   8,  10,
     *    9,   3,  18,  25,  53,  29,  43,   0,  17, 153,
     *   52,  55, 167,  44,  45, 182,  74,  46,  11, 187,
     *  201,  31,  27,  99, 209,  24,  22,  20,  21,  65,
     *   35,  19,  54,  23,   0,  14,  58,   0,  32,  12,
     *   47,  16,   0,  34,   0,   0,  30,   0,  13,  33,
     *   61,  63, 292,  57,  59,  66, 272,   0,  94, 175,
     *  226, 286,   0,   0,   0, 176, 227, 287,   0,   0,
     *    0,  96,   0, 177,   0, 267, 228, 288,   0,   0,
     *    0,   0,  93, 147, 162, 5*0,
     *    0,  50,   0,   0,   0,   0,  36,   0,  64,   0,
     *    0,  48,   0,   0, 148,   0,   0,  26,  49,  70,
     *  178,  97, 170, 229,   0, 180, 268, 230,   0, 289,
     *    0,   0,  15, 181,   0, 269, 4*0,
     *    0,   0,   0, 231,   0, 290,   0,  38,   0,   0,
     *  152,  39,  40,   0,  41, 232,   0, 291,   0,   0,
     *    0,   0,   0,  75, 154,   0,   0,   0, 183,   0,
     *    0,   0,   0,   0,   0,  98, 184, 234, 185, 270,
     *    0,   0,   0, 186,   0,   0, 271, 235,   0,   0,
     *   62,   0,   0,   0,   0,   0,   0, 101,   0, 188,
     *    0,   0,   0,   0,   0, 102, 189, 3*0,
     *  236,   0, 294,  67,   0, 190,   0,   0,   0, 295,
     *    0,   0, 104, 191, 237,   0, 105, 192, 274, 238,
     *  296, 112, 245, 303, 113, 199,   0, 278, 246,   0,
     *  304,   0,   0,   0,   0, 200,   0,   0, 279, 247,
     *    0, 305,   0,   0, 172, 5*0,
     *    0, 120, 122, 208,   0, 282, 255,   0, 312,   0,
     *  7*0,                               283, 256,   0,
     * 10*0,
     *  275, 194, 108, 241, 299, 202,   0,  68,  69,  71,
     *   72,  73,  42,  37,  76,  77,  78,  79,  80,  81,
     *   82,  83,  92,  95, 100, 103, 106, 107, 109, 110,
     *  111, 114, 115, 116, 117, 118, 119, 121, 123, 124,
     *  125, 126, 127, 128, 129, 149, 150, 151, 155, 156,
     *  157, 158, 159, 163, 164, 165, 166, 168, 169, 170,
     *  171, 193, 195, 196, 197, 198, 203, 204, 205, 206,
     *  207, 210, 211, 212, 213, 214, 215, 216, 217, 218,
     *  225, 233, 239, 240, 242, 243, 244, 248, 249, 250,
     *  251, 252, 253, 254, 257, 258, 259, 260, 262, 262,
     *  263, 264, 265, 266, 273, 276, 277, 280, 282, 284,
     *  285, 293, 297, 298, 300, 301, 302, 306, 307, 308,
     *  309, 310, 311,  60, 313, 314, 315, 316, 317, 318,
     *  319, 320, 321, 322, 323, 324,  84,  85,  86,  87,
     *   88,  89,  90,  91, 130, 131, 132, 133, 134, 135,
     *  136, 137, 138, 139, 140, 141, 142, 143, 144, 145,
     *  146, 160, 161, 173, 174, 210, 220, 221, 222, 223,
     *  224,16*0,  56/
c
      data iread /0/
c
c     call old Irwin routine MPARTF if desired
c
      if(irwtab.eq.0) then
         call mpartf(jatom,ion,indmol,t,u)
         return
      end if
c
c     read data if first call:
c
      if(iread.ne.1) then
        if(irwtab.eq.0) then
           open(67,file= './data/irwin_orig.dat',status='old')
         else
           open(67,file= './data/irwin_bc.dat',status='old')
        end if
        read(67,*)
        read(67,*)
        do j=1,92
          do i=1,3
            if(j.eq.1.and.i.eq.3) goto 10
            sp=float(j)+float(i-1)/100.
            read(67,*) spc,aa
            do k=1,6
                a(k,i,j)=aa(k)
            end do
   10    continue
         end do
       end do
c
       read(67,*)
       read(67,*)
       read(67,*)
       do i=1,324
          read(67,*,end=15) spec(i),aa
          do j=1,6
             am(j,i)=aa(j)
          end do
        end do
   15   continue
        close(67)
        iread=1
      endif
c
c     evaluation of the partition function
c     stop if T is out of limits of Irwin's tables
c
        if(t.lt.1000.) then
          stop 'partf; temp<1000 K'
        else if(t.gt.16000.) then
          stop 'partf; temp>16000 K'
        endif
        tl=log(t)
        u=0.
c
c     atomic species
c
      if(jatom.gt.0.and.ion.gt.0) then
        ulog=    a(1,ion,jatom)+
     *       tl*(a(2,ion,jatom)+
     *       tl*(a(3,ion,jatom)+
     *       tl*(a(4,ion,jatom)+
     *       tl*(a(5,ion,jatom)+
     *       tl*(a(6,ion,jatom))))))
       if(jatom.eq.5.and.ion.eq.3) ulog=1.
C      write(*,*) 'bor',ion,tl,ulog
c    * write(6,631) ion,tl,a(1,ion,jatom),tl*a(2,ion,jatom),
c    tl**2*a(3,ion,jatom),tl**3*a(4,ion,jatom),tl**4*a(5,ion,jatom),
c    * tl**5*a(6,ion,jatom),ulog
c 631  format('bor',i4,1p8e11.3)
        u=exp(ulog)
        return
      end if
c
c     molecular species
c
      if(indmol.gt.0) then
         indm=irwind(indmol)
         if(indm.le.0) return
         ulog=    am(1,indm)+
     *        tl*(am(2,indm)+
     *        tl*(am(3,indm)+
     *        tl*(am(4,indm)+
     *        tl*(am(5,indm)+
     *        tl*(am(6,indm))))))
        u=exp(ulog)
c       if(t.gt.5128..and.t.lt.5129.)
c    *  write(6,631) t,indmol,indm,u
c 631   format('irwpf',f10.1,2i5,f16.3)
      end if
      return
      end