科研技术基础

0. 前言

被 Tong Dang push学学Vim，记录一下学习的过程，争取能学下来（好难啊，什么阴间玩意）

update on 2023.2.6 ：经 Tong Dang 和 Binzheng Zhang 教导，似乎Vim不用学那么多，又自己高估对自己的要求了……

update on 2023.4.13 ：不想一遍遍翻 Tong Dang 的ppt了抽时间把超算的运行记录在本文里吧……

update on 2023.5.4 ：新增读取数据的代码（和我写的注释），生成视频

1. WSL

不想用辣鸡虚拟机，不想用学校的vlab，不想放弃windows……

那就用Windows Subsystem for Linux！

听说win11可以自带linux软件的GUI界面炒鸡香，果断更新。

（其实WSL似乎对俺这个大研的作用只有试试linux/vim和ffmpeg……）后来发现

需要记住的指令

1. 进入WSL（from Powershell）

   wsl

2. 查看WSL状态

   wsl -l -v

3. WSL关机

   wsl --shutdown

2. Linux/WSL

1. 从WSL返回到Powershell

   exit

2. 进入文件夹

   cd 'filename'

3. Vim

刚刚进入为 Normal 模式

1. 切换到 Insert 模式

   i

2. 回到 Normal 模式

   'Esc'

3. Normal 模式中，删除光标所在位置（光标线后）的字符

   x

4. 存盘

   :w 'filename'

5. 退出vim

   :q

6. 打开一个文件

   :e 'path'

7. 保存并退出

   :x

   ZZ

   :wq

8. 退出不保存

   :q!

4. 超算运行

咕咕咕

5. 数据读取

原代码来自 Tong Dang，俺给它改成木星的了，顺便加入注释供各位看官理解含义。

#python header file
import os
from collections import defaultdict
from time import time
import h5py
import numpy as np

from scipy.io import savemat


def get_equtorial_lfm(xCenter, yCenter, zCenter, value, kind='eq'):# Lyon-Fedder-Mobarry:the model of solar wind interacts with earth magsphere......
    
    xCenter, yCenter, zCenter, value = xCenter.T, yCenter.T, zCenter.T, value.T
    
    ni, nj, nk = zCenter.shape
    xe = np.zeros((ni, nj*2+1))
    ye = np.zeros_like(xe)
    ze = np.zeros_like(xe)
    data = np.zeros_like(xe)

    if kind == 'eq':# equitor 
        for i in range(ni):
            for j in range(nj):
                xe[i, j] = xCenter[i, j, 0]# z=0
                ye[i, j] = yCenter[i, j, 0]
                data[i, j] = (value[i, j, 0]+value[i, j, -1])/2

            for j in range(nj, nj*2):
                xe[i, j] = xCenter[i, 2*nj-1-j, nk//2-1]# z=180
                ye[i, j] = yCenter[i, 2*nj-1-j, nk//2-1]
                data[i, j] = (value[i, 2*nj-1-j, nk//2-1] + value[i, 2*nj-1-j, nk//2])/2

            xe[i, -1] = xe[i, 0]
            ye[i, -1] = ye[i, 0]
            data[i, -1] = data[i, 0]


    if kind == 'zx':# xz
        for i in range(ni):
            for j in range(nj):
                xe[i, j] = xCenter[i, j, nk//4]# z=90
                ye[i, j] = zCenter[i, j, nk//4]
                data[i, j] = value[i, j, nk//4]

            for j in range(nj, nj * 2):
                xe[i, j] = xCenter[i, 2 * nj-1 - j, nk // 4*3]# z=270
                ye[i, j] = zCenter[i, 2 * nj-1 - j, nk // 4*3]
                data[i, j] = value[i, 2 * nj-1 - j, nk // 4*3]

            xe[i, -1] = xe[i, 0]
            ye[i, -1] = ye[i, 0]
            data[i, -1] = data[i, 0]

    return xe, ye, ze, data

    pass

# file source
hdf_fldr = '/public/home/oybdooo/run_jupiter_0407/'# your h5 source
mat_dir  = '/public/home/oybdooo/matdata/jupiter_data_0412/'# the mat filepath u want to save

os.makedirs(mat_dir, exist_ok=True) # to create your mat data path (robustness)

i_mf = 1

file_list = os.listdir(hdf_fldr)
msphere_list = []
for c in file_list:
    if 'msphere' in c and 'Res' not in c:# to find the file we want to read eg:  msphere_0002_0002_0001_0000_0000_0000.gam.h5
        msphere_list.append(c)
msphere_list.sort()
msphere_list = msphere_list[1:]
prefix = msphere_list[-1][:-21]# get msphere_0002_0002_0001_0000 to msphere[-1](abandon the last 21 char)
name_split = msphere_list[-1].split('_')# divide prefix with _
NI_rank = int(name_split[1])# 0002
NJ_rank = int(name_split[2])# 0002
NK_rank = int(name_split[3])# 0001
N_rank = int(name_split[-1].split('.')[0])# 0000

# ------------------------------------------------------------------
if i_mf == 1:# judge whether it is multi-fluid
    i_mf_list = ['D1', 'D2', 'P1', 'P2', 'Vx1', 'Vx2', 'Vy1', 'Vy2',
                 'Vz1', 'Vz2']
else:
    i_mf_list = []

norm_list = ['D', 'Vx', 'Vy', 'Vz', 'P', 'Bx', 'By', 'Bz', 'Jx', 'Jy', 'Jz']+i_mf_list

first_list = ['X', 'Y', 'Z', 'dV','X1']# dV: volume element
# ------------------------------------------------------------------


Data_merge_first = defaultdict(np.ndarray)

step_start = 0
step_end   = 2400 # run for every step
for step in range(step_start, step_end,1):
    t0 = time()
    Data_merge = defaultdict(np.ndarray)
    Data_lon = defaultdict(list)
    rank_count=0
    for irank in range(0,NI_rank):
        Data = defaultdict(list)
        for jrank in range(0,NJ_rank):
            fn = prefix+'{:0>4d}_{:0>4d}_0000.gam.h5'.format(irank,jrank)# find the file
            rank_count = rank_count+1
            file_up = h5py.File(hdf_fldr+'/'+fn, 'r')# read the file
            file = file_up['Step#{}'.format(step)]
            if step == step_start:# only gird file need save x,y,z!
                x = np.array(file_up['X'], order='C')
                y = np.array(file_up['Y'], order='C')
                z = np.array(file_up['Z'], order='C')
                v = np.array(file_up['dV'], order='C')
                nip, njp, nkp = x.shape
                I = np.arange(1, nip)# [1,2,3,...,nip]
                J = np.arange(1, njp)
                K = np.arange(1, nkp)

                def center_3d(v):# find the center point
                    return 0.125*(v[:-1, :-1, :-1]+v[1:, :-1, :-1]+v[:-1, 1:, :-1]+v[:-1, :-1, 1:]
                                  +v[1:, 1:, :-1]+v[1:, :-1, 1:]+v[:-1, 1:, 1:]+v[1:, 1:, 1:])

                xc = center_3d(x)
                yc = center_3d(y)
                zc = center_3d(z)

                Data['X'].append(xc)
                Data['Y'].append(yc)
                Data['Z'].append(zc)
                Data['dV'].append(v)
                Data['X1'].append(x)


            for c in norm_list:
                Data[c].append(np.array(file[c], order='C'))
        if step == step_start:# only gird file need save x,y,z!
            for c in first_list:
                Data_lon[c].append(np.hstack(Data[c]))

        for c in norm_list:
            Data_lon[c].append(np.hstack(Data[c]))

    if step == step_start:
        for c in first_list:
            Data_merge_first[c] = np.dstack(Data_lon[c])

    for c in norm_list:
        Data_merge[c] = np.dstack(Data_lon[c])

    def getdata(kind='eq'):
        Data_interp = defaultdict(np.ndarray)
        for c in norm_list:
            _, _, _, res = get_equtorial_lfm(Data_merge_first['X'], Data_merge_first['Y'], Data_merge_first['Z'], Data_merge[c], kind)
            Data_interp[c] = res
        return Data_interp

    def name_transform(name, label='e'):# transform name rhoe/rhom blablabla
        name = name.lower().replace('d', 'rho')
        if name[-1].isdigit():
            name = name[:-1]+label+name[-1]
        else:
            name = name+label
        return name

    # get and save equatorial data
    Data_equator = getdata('eq')
    save_path_equator = mat_dir+'/eq-{:0>4d}.mat'.format(step)
    save_dict_equator = defaultdict(np.ndarray)
    for c in norm_list:
        save_dict_equator[name_transform(c, label='e')] = Data_equator[c]
    xe, ye, ze, _ = get_equtorial_lfm(Data_merge_first['X'], Data_merge_first['Y'], Data_merge_first['Z'], Data_merge['D'], 'eq')
    savemat(save_path_equator, save_dict_equator, do_compression=True)

    # get and save xz data
    Data_xz = getdata('zx')
    save_path_xz = mat_dir+'/xz-{:0>4d}.mat'.format(step)
    save_dict_xz = defaultdict(np.ndarray)
    for c in norm_list:
        save_dict_xz[name_transform(c, label='m')] = Data_xz[c]
    xm, ym, zm, _ = get_equtorial_lfm(Data_merge_first['X'], Data_merge_first['Y'], Data_merge_first['Z'], Data_merge['D'], 'zx')
    savemat(save_path_xz, save_dict_xz, do_compression=True)
    print('Step: {} time usage: {}'.format(step, time()-t0))    
    
    if step == step_start:
        # save grid
        savemat(mat_dir+'/grid_eq.mat', dict(xe=xe, ye=ye, ze=ze), do_compression=True)
        savemat(mat_dir+'/grid_xz.mat', dict(xm=xm, ym=ym, zm=zm), do_compression=True)

6. 生成视频

俺是先用WSL进入图片对应路径，如果图片的长度和宽度像素不为偶数需要用scale调整到偶数！

ffmpeg -framerate 25 -pattern_type glob -i '*.png' -vf scale=2364:2262 -pix_fmt yuv420p eq-Bt.mp4