物理の駅 Physics station by 現役研究者

テクノロジーは共有されてこそ栄える

LISE++ファイルをPythonで読み取りBeamDumpの設定を調べる

LISE++ファイルは、テキストファイルで作られているのでPythonで容易に読み取ることができる。 [] で囲まれた文字列をカテゴリ、その中で記述された = で区切られたKeyとValueを辞書型で詰めている。; 以降はコメントとして無視している。

def read_lise_lpp(filename):
    import re
    try:
        f = open(filename)
        lines = f.readlines()
    except:
        f = open(filename,encoding="utf-8")
        lines = f.readlines()
    f.close()

    obj = {}
    obj["Version"] = lines[0].strip()
    catgegory = None
    for line in lines[1:]:
        line = line.strip()
        if line=="":continue

        # section
        if line[0]=="{" and line[-1]=="}":continue
        
        # catgegory
        if line[0]=="[" and line[-1]=="]":
            catgegory=line[1:-1]
            obj[catgegory]={}
            continue
        
        # comments
        if ';' in line:
            comment_pos = line.index(';')
            line = line[:comment_pos]
            
        # key and value
        line = line.strip(' ')
        equal_pos = line.index('=')
        key = line[0:equal_pos]
        value = line[equal_pos+1:]
        key = key.strip(' ')
        key = re.sub('\s+', ' ', key)
        value = value.strip(" ")
        value = re.sub('\s+', ' ', value)
        obj[catgegory][key]=value
    return obj

RIBFのBigRIPSで実験するときのBeamDumpの設定を調べる関数を紹介する。RIBFのウェブサイトにPOSTでデータ取得している。

def get_beam_dump(obj,verbose=False):
    AZQ = obj["settings"]["A,Z,Q"]
    Energy = obj["settings"]["Energy"]
    TargetContents = obj["target"]["Target contents"]
    TargetThickness = obj["target"]["Target thickness"].split(",")[1]
    D1Brho = float(obj["D1_DipoleSettings"]["Brho"].replace(" Tm",""))
    return get_beam_dump_impl(AZQ,Energy,TargetContents,TargetThickness,D1Brho,verbose)
    
def get_beam_dump_impl(AZQ,Energy,TargetContents,TargetThickness,D1Brho,verbose=False):
    import requests
    import json
    arguments={"Status":"Error","AZQ":AZQ,
               "Energy":Energy,
               "TargetContents":TargetContents,
               "TargetThickness":TargetThickness,
               "D1Brho":D1Brho}
    if   AZQ=="238U 86+" and Energy=="345 MeV/u":Bnum=0
    elif AZQ=="136Xe52+" and Energy=="345 MeV/u":Bnum=1
    elif AZQ=="124Xe52+" and Energy=="345 MeV/u":Bnum=2
    elif AZQ=="86Kr36+" and Energy=="345 MeV/u":Bnum=3
    elif AZQ=="78Kr36+" and Energy=="345 MeV/u":Bnum=4
    elif AZQ=="70Zn30+" and Energy=="345 MeV/u":Bnum=5
    elif AZQ=="48Ca20+" and Energy=="345 MeV/u":Bnum=6
    elif AZQ=="40Ca20+" and Energy=="345 MeV/u":Bnum=7
    elif AZQ=="36Ar18+" and Energy=="345 MeV/u":Bnum=8
    elif AZQ=="18O 8+" and Energy=="345 MeV/u":Bnum=9
    elif AZQ=="18O 8+" and Energy=="250 MeV/u":Bnum=10
    elif AZQ=="18O 8+" and Energy=="230 MeV/u":Bnum=11
    elif AZQ=="16O 8+" and Energy=="345 MeV/u":Bnum=12
    elif AZQ=="16O 8+" and Energy=="250 MeV/u":Bnum=13
    elif AZQ=="14N 7+" and Energy=="250 MeV/u":Bnum=14
    elif AZQ=="4He2+" and Energy=="320 MeV/u":Bnum=15
    elif AZQ=="2H 1+" and Energy=="250 MeV/u":Bnum=16
    else:
        if verbose:print(arguments)
        return arguments
    if type(TargetThickness) is not str:
        TargetThickness = str(TargetThickness)
    if TargetContents=="0,4,1,9.012" or TargetContents=="Be":
        if TargetThickness=="0.1":Tnum=0
        elif TargetThickness=="0.5":Tnum=1
        elif TargetThickness in ["1","2","3","4","5","6","7","8","9","10","11","12"]:Tnum=int(TargetThickness)+2
        elif TargetThickness=="15":Tnum=15
        elif TargetThickness=="17":Tnum=16
        elif TargetThickness=="20":Tnum=17
        elif TargetThickness=="30":Tnum=18
        elif TargetThickness=="40":Tnum=19
        else:
            if verbose:print(arguments)
            return arguments
    elif TargetContents=="0,74,1,183.85" or TargetContents=="W":
        if TargetThickness=="0.1":Tnum=20
        elif TargetThickness=="0.2":Tnum=21
        elif TargetThickness=="0.3":Tnum=22
        elif TargetThickness=="0.5":Tnum=23
        elif TargetThickness=="1":Tnum=24
        elif TargetThickness=="1.4":Tnum=25
        elif TargetThickness=="2.1":Tnum=26
        elif TargetThickness=="2.8":Tnum=27
        else:
            if verbose:print(arguments)
            return arguments
    else:
        if verbose:print(arguments)
        return arguments
    
    if TargetContents=='0,4,1,9.012':
         TargetContents='Be'
    elif TargetContents=='0,74,1,183.85':
        TargetContents='W'
    payload = {'Bnum': Bnum,'Tnum': Tnum,'D1Brho': D1Brho,'btnExec': ''}
    files = {(None, None)}
    ret = requests.post('https://ribf.riken.jp/BigRIPSInfo/beamdump/result.php', files=files,data=payload)
    lines = ret.text.split("\n")
    status={}
    status["Status"]="OK"
    left = lines[46][lines[46].index(">")+1:lines[46].index('mm')].strip(' ')
    right = lines[48][lines[48].index(">")+1:lines[48].index('mm')].strip(' ')
    status["Left"]= "---" if left=="---" else float(left)
    status["Right"]= "---" if right=="---" else float(right)
    if verbose:
        payload['Beam']=AZQ+' '+Energy
        payload['Target']=TargetContents+' '+TargetThickness+'mmt'
        del payload['btnExec']
        print(payload,"-->",status)
    return status

使い方は直接パラメータを指定するか、LISE++ファイルから読み取るかの2種類ある。

# 直接パラメータを指定する方法
print(get_beam_dump_impl(AZQ="238U 86+",Energy="345 MeV/u",TargetContents="Be",TargetThickness=0.1,D1Brho=7.0,verbose=True))

# LISEファイルから関数`read_lise_lpp`で読み取った辞書型を使う方法
import requests
url = "http://ribf.riken.jp/BigRIPSInfo/lise/238U_Be_fission_Standard_LAA.lpp"
with open("lise.lpp", "w") as f:
    f.write(requests.get(url).text)
print(get_beam_dump(read_lise_lpp("lise.lpp"),verbose=True))

標準出力は次のようになる

{'Bnum': 0, 'Tnum': 0, 'D1Brho': 7.0, 'Beam': '238U 86+ 345 MeV/u', 'Target': 'Be 0.1mmt'} --> {'Status': 'OK', 'Left': -55.0, 'Right': 125.0}
{'Status': 'OK', 'Left': -55.0, 'Right': 125.0}

{'Bnum': 0, 'Tnum': 6, 'D1Brho': 8.1005, 'Beam': '238U 86+ 345 MeV/u', 'Target': 'Be 4mmt'} --> {'Status': 'OK', 'Left': -125.0, 'Right': 125.0}
{'Status': 'OK', 'Left': -125.0, 'Right': 125.0}

RIBFの Calculation for the position of the beam dump で直接調べても良い (が大量にチェックする時は大変)。まあ、用途は限られているが...

Python jupyter lab で動的なグラフを作る

Jupyter lab等でグラフをインタラクティブに動かしたいとき、例えば特定の値を変えたグラフをセルを実行することなく次々と表示させたい時、ipywidgets.interactを使うと良い。

使用例

%matplotlib inline
import matplotlib.pyplot as plt
import ipywidgets
import numpy as np

plt.rcParams['figure.figsize'] = [5, 4]

def f(seed):
    np.random.seed(seed)
    plt.hist(np.random.randn(1000),range=[-5,5],bins=50)
    plt.ylim(0,110)
    plt.show()

ipywidgets.interact(f, seed=(0,100,1))

seedを0から100まで変えたときの、正規分布を確認できるだろう。

もっと詳しく

qiita.com

Python 時間をdatetime型、文字列、unixtimeに相互変換

文字列→datetime, unixtime, 文字列

import datetime

tstr = '2020/9/17 9:00:01'
tdt = datetime.datetime.strptime(tstr, '%Y/%m/%d %H:%M:%S')
d_week = {'Sun':'日','Mon':'月','Tue':'火','Wed':'水','Thu':'木','Fri':'金','Sat':'土'}
print(type(tdt))
print(tdt)
print(tdt.timestamp())
print(tdt.strftime('%Y年%m月%d日({}) %H時%M分%S秒').format(d_week[tdt.strftime('%a')]))
print(tdt.strftime('%Y/%m/%d %H:%M:%S'))

unixtime→datetime, unixtime, 文字列

import datetime

tdt = datetime.datetime.fromtimestamp(1600300801.0)
d_week = {'Sun':'日','Mon':'月','Tue':'火','Wed':'水','Thu':'木','Fri':'金','Sat':'土'}
print(type(tdt))
print(tdt)
print(datetime.datetime.timestamp(tdt))
print(tdt.strftime('%Y年%m月%d日({}) %H時%M分%S秒').format(d_week[tdt.strftime('%a')]))
print(tdt.strftime('%Y/%m/%d %H:%M:%S'))

datetime→datetime, unixtime, 文字列

import datetime

tdt = datetime.datetime(2022, 9, 17, 9, 0, 1)
d_week = {'Sun':'日','Mon':'月','Tue':'火','Wed':'水','Thu':'木','Fri':'金','Sat':'土'}
print(type(tdt))
print(tdt)
print(datetime.datetime.timestamp(tdt))
print(tdt.strftime('%Y年%m月%d日({}) %H時%M分%S秒').format(d_week[tdt.strftime('%a')]))
print(tdt.strftime('%Y/%m/%d %H:%M:%S'))

出力

<class 'datetime.datetime'>
2022-09-17 09:00:01
1663372801.0
2022年09月17日(土) 09時00分01秒
2022/09/17 09:00:01

現在時刻から

tdt = datetime.datetime.now()

matplotlibで時間のフォーマットを変える

import matplotlib.dates as mdates
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%d日%H時'))

ツイッター用

dt = datetime.datetime.strptime(data["created_at"],'%a %b %d %H:%M:%S %z %Y')
JST = datetime.timezone(datetime.timedelta(hours=9), "JST")
dt = dt.astimezone(JST)
diff_days = (datetime.datetime.now(JST)-dt).days

import datetime

#tweet IDをタイムスタンプに
def tweet_id2timestamp(tweet_id):
    return ((tweet_id >> 22) + 1288834974657) / 1000

#tweet IDをdatetimeに
def tweet_id2time(tweet_id):
    return datetime.datetime.fromtimestamp(tweet_id2timestamp(tweet_id))

#datetimeをtweet_idに
def time2tweet_id(dt):
    return timestamp2tweet_id(dt.timestamp())

#timestampをtweet_idに
def timestamp2tweet_id(timestamp):
    return (round(timestamp*1000)-1288834974657)<<22

その日が月末かどうかを判定する方法

import datetime

def is_end_of_month(date):
    first_day_of_month = date + datetime.timedelta(days=1)
    return date.month != first_day_of_month.month

date1 = datetime.datetime(2021, 12, 31) #年末
date2 = datetime.datetime(2020, 2, 29) #閏年
date3 = datetime.datetime(2022, 3, 15) #普通の日

print(is_end_of_month(date1))
print(is_end_of_month(date2))
print(is_end_of_month(date3))

検索ワード: 曜日

Python ログインが必要なサイトでスクレイピングするための前準備 POSTで認証

session を使い回すところがポイント

import requests

session = requests.Session()

url = "https://test.test/login"
data = {
    "email": "mail@mail.mail",
    "password": "password",
}

ret = session.post(url=url, params=data)
assert(ret.status_code==200)

ret = session.get("https://test.test/data")
with open("data.data", "w") as f:
    f.write(ret.text)

(Solved v2.2.2 and later) Slackdump_v2.2.0 Katakana charactors directory problem

日本語訳は下に (Japanese translation below)

This problem was fixed in version v2.2.2 and later. Please download the latest version of the slackdump.

There is an excelent software called slackdump that allows you to back up all your visible Slack messages and attachments without special permissions, i.e. administrator or owner.

github.com

When the channel name contains Japanese full-width katakana characters, the directory name to save the slack message data (eg. 2022-10-10.json) will be half-width katakana with using slackdump. Attachments are saved in correct full-width katakana directories.

To temporarily solve this problem, I prepared a python script that copies the message data = JSON files in the half-width katakana directory to the full-width katakana directory.

import zipfile
import mojimoji # need install mojimoji to convert half-width katakana to full-width katakana
import json
import shutil
import os
zipPath_org = r"C:\Users\Masahiro\Downloads\escan.zip"
zipPath = os.path.join(os.path.dirname(zipPath_org),
                       os.path.splitext(os.path.basename(zipPath_org))[0]+"_new"+os.path.splitext(os.path.basename(zipPath_org))[1])

shutil.copy(zipPath_org,zipPath)
with zipfile.ZipFile(zipPath, 'a') as myzip:
    infolist = myzip.infolist()

    channels = json.load(myzip.open('channels.json'))
    channel_names = []
    for channel in channels:
        channel_names.append(channel['name'])
        # print(channel['name']) for debug
        
    for info in infolist:
        if info.filename.split('/')[0] not in channel_names:
            channel_half_width = info.filename.split('/')[0]
            channel_full_width = mojimoji.han_to_zen(channel_half_width, ascii=False, kana=True, digit=False)
            if channel_half_width == channel_full_width:continue
            if channel_full_width not in channel_names:
                print(channel_full_width)
            assert(channel_full_width in channel_names)
            filename_new = channel_full_width+info.filename[len(channel_half_width):]
            # print(info.filename,"->",filename_new) # for debug
            myzip.writestr(filename_new,myzip.open(info.filename).read())
    myzip.close()
print(zipPath, "is output.")

(Japanese translation)

バージョンv2.2.2以降で、この問題は修正済みです最新バージョンのslackdumpをダウンロードしてください。

Slack上で見えているメッセージを、管理者・オーナーなどの特別な権限がなくても、バックアップできる slackdumpという素晴らしいソフトウェアを見つけました。

github.com

チャンネル名に日本語の全角カタカナ文字を含む場合、slackdumpでメッセージデータ(例: 2022-10-10.json)の保存先のディレクトリ名が、半角カタカナになってしまうバグがあるようです。添付ファイルは、全角カタカナのディレクトリに正しく保存されています。

この問題を一時的に解決するため、半角カタカナのディレクトリにあるメッセージデータ = JSONファイルを全角カタカナ文字のディレクトリにコピーするPythonスクリプトを用意しました。

参照

kenkyu-note.hatenablog.com