Upload New File

import numpy as np

import h5py

def load_data(filedir, n0, n1):

    # Set minimum value for the magnetic field and velocity

    epsilon = 1e-6

    # Open HDF5 file

    filename = filedir + 'full_dtb_R_007'

    infile = h5py.File(filename,'r')

    # Define datasets

    dsetrho = "Density"

    dsetT = "Temperature"

    dsetvx = "x-velocity"

    dsetvy = "y-velocity"

    dsetvz = "z-velocity"

    dsetBx = "Bx"

    dsetBy = "By"

    dsetBz = "Bz"

    # Define cutout corner and size

    nx0 = n0[0]

    ny0 = n0[1]

    nz0 = n0[2]

    nx = n1[0]

    ny = n1[1]

    nz = n1[2]

    # Read input data, subset and normalise

    rhomin = -2.8

    rhomax = 2.5

    dsetaux = infile[dsetrho]

    rhoaux = np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz])

    rho = np.log10(rhoaux)

    #rhomin = np.amin(rho)

    #rhomax = np.amax(rho)

    #print("Density: ",rhomin,rhomax)

    rhoaux = (rho - rhomin)/(rhomax-rhomin)

    rho = np.clip(rhoaux,0.0,1.0)

    rhomin = np.amin(rho)

    rhomax = np.amax(rho)

    print("Density: ",rhomin,rhomax)

    rhomin = 0.0 

    rhomax = 9.0

    dsetaux = infile[dsetT]

    Taux = np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz])

    T = np.log10(Taux)

    #print("Temperature: ",rhomin,rhomax)

    Taux = (T - rhomin)/(rhomax-rhomin)

    T = np.clip(Taux,0.0,1.0)

    rhomin = np.amin(T)

    rhomax = np.amax(T)

    print("Temperature: ",rhomin,rhomax)

    rhomin = np.log10(epsilon)

    rhomax = np.log10(0.1)

    dsetaux = infile[dsetBx]

    B2 = np.square(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))

    Bx = np.log10(np.abs(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))+epsilon)

    dsetaux = infile[dsetBy]

    B2 = B2 + np.square(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))

    By = np.log10(np.abs(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))+epsilon)

    dsetaux = infile[dsetBz]

    B2 = B2 + np.square(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))

    Bz = np.log10(np.abs(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))+epsilon)

    B2 = np.log10(np.sqrt(B2)+epsilon)

    Baux = (B2 - rhomin)/(rhomax-rhomin)

    B2 = np.clip(Baux,0.0,1.0)

    Baux = (Bx - rhomin)/(rhomax-rhomin)

    Bx = np.clip(Baux,0.0,1.0)

    Baux = (By - rhomin)/(rhomax-rhomin)

    By = np.clip(Baux,0.0,1.0)

    Baux = (Bz - rhomin)/(rhomax-rhomin)

    Bz = np.clip(Baux,0.0,1.0)

    print("B: ",np.amin(B2),np.amin(Bx),np.amin(By),np.amin(Bz))

    print("B: ",np.amax(B2),np.amax(Bx),np.amax(By),np.amax(Bz))

    infile.close()

    filename = filedir + 'full_v_R_007'

    infile = h5py.File(filename,'r')

    rhomin = 0.0

    rhomax = 0.05

    dsetaux = infile[dsetvx]

    #print("vx: ",np.amin(dsetaux),np.amax(dsetaux))

    vx = np.abs(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))

    dsetaux = infile[dsetvy]

    #print("vy: ",np.amin(dsetaux),np.amax(dsetaux))

    vy = np.abs(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))

    dsetaux = infile[dsetvz]

    #print("vz: ",np.amin(dsetaux),np.amax(dsetaux))

    vz = np.abs(np.array(dsetaux[nx0:nx,ny0:ny,nz0:nz]))

    vaux = (vx - rhomin)/(rhomax-rhomin)

    vx = np.clip(vaux,0.0,1.0)

    vaux = (vy - rhomin)/(rhomax-rhomin)

    vy = np.clip(vaux,0.0,1.0)

    vaux = (vz - rhomin)/(rhomax-rhomin)

    vz = np.clip(vaux,0.0,1.0)

    print("velocity: ",np.amin(vx),np.amin(vy),np.amin(vz))

    print("velocity: ",np.amax(vx),np.amax(vy),np.amax(vz))

    infile.close()

    # Output data for debugging

    #Bz.tofile('data.bin')

    return [rho, T, vx, vy, vz, Bx, By, Bz, B2]