Environment

GeForce GTX 1070 (8GB)
ASRock Z170M Pro4S [Intel Z170chipset]
Ubuntu 16.04 LTS desktop amd64
TensorFlow v1.2.1
cuDNN v5.1 for Linux
CUDA v8.0
Python 3.5.2
IPython 6.0.0 -- An enhanced Interactive Python.
gcc (Ubuntu 5.4.0-6ubuntu1~16.04.4) 5.4.0 20160609
GNU bash, version 4.3.48(1)-release (x86_64-pc-linux-gnu)
scipy v0.19.1
geopandas v0.3.0
MATLAB R2017b (Home Edition)
ADDA v.1.3b6

Bug in previous versions (up to v0.3)

[obsolete] bisphere.f > tool > calc_qext_bisphere_180104.py > v0.1-v0.3 > calculate Qext from the result (R1, R2, CEXT)

In the v0.3, the code reads [result.txt] not [bisphere.print] mistakenly.

v0.4

read [bisphere.print]
- instead of [result.txt]

Format difference

There is a difference in the format between [result.txt] and [bisphere.txt] in terms of the expression of the values as shown below (e.g. R1).

result.txt(excerpt)

This test result was computed on an IBM RISC model 397 workstation
for the current code setting.

R1= .5000D+01  R2= .5000D+01  R12= .2000D+02  LAM= .6283D+01
X12= .2000D+02  NODRT= 33
X1= .500000D+01 N1= .150000D+01 K1= .500000D-02 NODR(1)= 14
X2= .500000D+01 N2= .150000D+01 K2= .500000D-02 NODR(2)= 14
 CEXT= .589017D+03  CSCA= .567076D+03    W= .962749D+00  <COS>= .717082D+00
 TEST OF VAN DER MEE & HOVENIER IS SATISFIED

bisphere.print(excerpt)

R1=0.5000D+01  R2=0.5000D+01  R12=0.1000D+02  LAM=0.6283D+01
X12=0.1000D+02  NODRT= 21
X1=0.500000D+01 N1=0.150000D+01 K1=0.500000D-02 NODR(1)= 14
X2=0.500000D+01 N2=0.150000D+01 K2=0.500000D-02 NODR(2)= 14
 CEXT=0.516661D+03  CSCA=0.495226D+03    W=0.958513D+00  <COS>=0.728258D+00
 TEST OF VAN DER MEE & HOVENIER IS SATISFIED
 
  S      ALPHA1      ALPHA2      ALPHA3      ALPHA4       BETA1       BETA2
  0     1.00000     0.00000     0.00000     0.93639     0.00000     0.00000
  1     2.18478     0.00000     0.00000     2.23908     0.00000     0.00000

In v0.4, the code reads the format for [bisphere.print].

code v0.4

calc_qext_bisphere_180104.py

import numpy as np
import sys

'''
v0.4 Jan. 06, 2018
  - read_param_fromfile() splits with delimiters ('=' and ' ')
  - fix bug: read [bisphere.print]
    + not [result.txt]
v0.3 Jan. 04, 2018
  - fix bug: calc_qext() uses reff not reff^2
v0.2 Jan. 04, 2018
  - rename to [calc_qext_bisphere_180104.py]
    + was [calc_qext_180104.py]
  - read_param_fromfile() uses convert_to_value()
  - add convert_to_value()
  - add calc_qext()
  - add calc_reff()
v0.1 Jan. 04, 2018
  - add read_param_fromfile()
'''


def convert_to_value(valstr):
    # to avoid:
    # TypeError: unsupported operand type(s) for
    # ** or pow(): 'str' and 'int'
    work = valstr.replace("D", "E")
    return np.float64(work)


def read_param_fromfile(akey):
    # 1. find line
    keyline = []
    with open('bisphere.print') as fd:
        lines = fd.readlines()
        for aline in lines:
            if not akey + "=" in aline:
                continue
            keyline += [aline]
    if len(keyline) == 0:
        return None
    # 2. extract value next to the keyword
    found = False
    txt = keyline[0].replace('=', ' ')
    for elem in txt.split(' '):
        if found:
            return convert_to_value(elem)
        found = akey in elem
    return None  # just in case


def calc_reff(r1, r2):
    # volume equivalent sphere radius
    return (r1**3 + r2**3)**(1./3.)


def calc_qext(r1, r2, cext):
    reff = calc_reff(r1, r2)
    denom = np.pi * reff * reff
    if abs(denom < sys.float_info.epsilon):
        return 0.0
    return cext / denom


r1 = read_param_fromfile("R1")
r2 = read_param_fromfile("R2")
cext = read_param_fromfile("CEXT")
reff = calc_reff(r1, r2)

qext = calc_qext(r1, r2, cext)
print("r1, r2, reff, cext")
print("%.3e %.3e %.5e %.3e" % (r1, r2, reff, cext))
print("qext")
print("%.7e" % qext)

run

$ python3 calc_qext_bisphere_180104.py 
r1, r2, reff, cext
5.000e+00 5.000e+00 6.29961e+00 5.167e+02
qext
4.1440896e+00