How to use the spacepy.datamodel.dmarray function in spacepy

Alpha = [alpha]

    # required setup
    MagEphemInfo = Lgm_MagEphemInfo.Lgm_MagEphemInfo(len(Alpha), 0)

    # setup a shortcut to MagModelInfo
    mmi = MagEphemInfo.LstarInfo.contents.mInfo.contents
    Lgm_Set_Coord_Transforms( datelong, utc, mmi.c) # dont need pointer as it is one

    # convert to **GSM**
    if coord_system == 'GSM':
        try:
            Pgsm = Lgm_Vector.Lgm_Vector(*pos)
        except TypeError:
            raise(TypeError("Position must be listlike" ) )
        ans['position']['GSM'] = datamodel.dmarray(pos, attrs={'units':'Re'})
    elif coord_system == 'SM':
        try:
            Psm = Lgm_Vector.Lgm_Vector(*pos)
        except TypeError:
            raise(TypeError("Position must be listlike" ) )
        Pgsm = Lgm_Vector.Lgm_Vector()
        Lgm_Convert_Coords( pointer(Psm), pointer(Pgsm), SM_TO_GSM, mmi.c )
        ans['position']['SM'] = datamodel.dmarray(pos, attrs={'units':'Re'})
        ans['position']['GSM'] = datamodel.dmarray(Pgsm.tolist(), attrs={'units':'Re'})
    else:
        raise(NotImplementedError("Only GSM or SM input currently supported"))

    Pwgs = Lgm_Vector.Lgm_Vector()
    Pmlt = Lgm_Vector.Lgm_Vector()
    Lgm_Convert_Coords( pointer(Pgsm), pointer(Pwgs), GSM_TO_WGS84, mmi.c )
    Lgm_Convert_Coords( pointer(Pwgs), pointer(Pmlt), WGS84_TO_EDMAG, mmi.c )

def runPA(quality):
    #test for range of pitch angles
    ticks = spt.tickrange('2002-04-18', '2002-04-19', 1)
    loci = spc.Coords([[-4,0,0], [-5,0,0]], 'SM', 'car')
    vals = []
    for pp in range(1,91):
        dum = lgmpy.get_Lstar(loci.data[1], ticks.UTC[1], alpha=pp, coord_system='SM', 
                              Bfield='Lgm_B_cdip', extended_out=True, LstarQuality=quality)
        try:
            vals.extend(dum[pp]['Lstar'])
        except (IndexError, TypeError): #get_Lstar returns a 0-D nan in some cases...
            vals.extend([dum[pp]['Lstar'].tolist()])
    fig = plt.figure()
    expect = float(tb.hypot(loci.data[1]))
    plt.plot(100*(dm.dmarray(vals)-expect)/expect, drawstyle='steps-mid')
    plt.ylim([-0.01, 0.01])
    plt.xlabel('Pitch Angle')
    plt.ylabel('100*(L* - L*$_{exp}$)/L*$_{exp}$')
    plt.title('Cent. dipole [-5,0,0]$_{SM}$'+'; Quality ({0})'.format(str(quality)))
    figname = 'MagEphem_CDIP_test_q{0}'.format(str(quality))   #'LGM_L*_vs_PA_cdip_q{0}'.format(str(quality))
    plt.savefig(''.join([figname,'.png']), dpi=300)
    plt.savefig(''.join([figname,'.pdf']))

raise(TypeError("Position must be listlike" ) )
        Pgsm = Lgm_Vector.Lgm_Vector()
        Lgm_Convert_Coords( pointer(Psm), pointer(Pgsm), SM_TO_GSM, mmi.c )
        ans['position']['SM'] = datamodel.dmarray(pos, attrs={'units':'Re'})
        ans['position']['GSM'] = datamodel.dmarray(Pgsm.tolist(), attrs={'units':'Re'})
    else:
        raise(NotImplementedError("Only GSM or SM input currently supported"))

    Pwgs = Lgm_Vector.Lgm_Vector()
    Pmlt = Lgm_Vector.Lgm_Vector()
    Lgm_Convert_Coords( pointer(Pgsm), pointer(Pwgs), GSM_TO_WGS84, mmi.c )
    Lgm_Convert_Coords( pointer(Pwgs), pointer(Pmlt), WGS84_TO_EDMAG, mmi.c )
    R, MLat, MLon, MLT = c_double(), c_double(), c_double(), c_double(),
    Lgm_EDMAG_to_R_MLAT_MLON_MLT( pointer(Pmlt),  pointer(R), pointer(MLat), pointer(MLon),
        pointer(MLT), mmi.c)
    ans['MLT'] = datamodel.dmarray(MLT.value, attrs={'coord_system': 'EDMAG'})

    # save Kp
    # TODO maybe add some Kp checking
    ans['Kp'] = datamodel.dmarray([Kp])

    # Set the LstarQuality, TODO add a comment here on what each does
    MagEphemInfo.LstarQuality = LstarQuality

    # L* in one place is L* in lots of places (for GPS set to False)
    MagEphemInfo.SaveShellLines = extended_out
    # TODO maybe not hardcoded, but for now its fine
    MagEphemInfo.LstarInfo.contents.VerbosityLevel = 0
    MagEphemInfo.LstarInfo.contents.mInfo.contents.VerbosityLevel = 0

    #MagEphemInfo->LstarInfo->mInfo->Bfield        = Lgm_B_T89;
    #MagEphemInfo->LstarInfo->mInfo->Bfield        = Lgm_B_cdip;

# Save nAlpha, and Alpha array to MagEphemInfo structure
    MagEphemInfo.nAlpha = len(Alpha)
    for i in range(len(Alpha)):
        MagEphemInfo.Alpha[i] = Alpha[i]

    # Set Tolerances
    Lgm_SetLstarTolerances(LstarQuality, 24, MagEphemInfo.LstarInfo )
    # *  Blocal at sat location
    MagEphemInfo.P = Pgsm

    Bvec = Lgm_Vector.Lgm_Vector()
    # Get B at the point in question
    MagEphemInfo.LstarInfo.contents.mInfo.contents.Bfield(pointer(Pgsm), pointer(Bvec),
                                                          MagEphemInfo.LstarInfo.contents.mInfo)
    ans['Bcalc'] = datamodel.dmarray(Bvec.tolist(), attrs={'units':'nT'})
    ans['Bcalc'].attrs['model'] = Bfield
    ans['Bcalc'].attrs['params'] = params
    ans['Bcalc'].attrs['coord_system'] = 'GSM'

    # save its magnitude in the structure
    MagEphemInfo.B = Bvec.magnitude()

    # check and see if the field line is closed before doing much work
    trace, northern, southern, minB, Lsimple = Closed_Field.Closed_Field(MagEphemInfo, extended_out=True)

    # presetup the ans[Angle] so that it can be filled correctly
    for pa in Alpha:
        ans[pa] = datamodel.SpaceData()
        ans[pa]['Lsimple'] = datamodel.dmarray([Lsimple])
        #sets up nans in case of Lstar failure
        ans[pa]['I'] = datamodel.dmarray(numpy.nan)

ans[pa]['Lstar'] = datamodel.dmarray(numpy.nan, attrs={'info':trace})
        ans[pa]['LMcIlwain'] = datamodel.dmarray(numpy.nan)
        ans[pa]['LHilton'] = datamodel.dmarray(numpy.nan)
        ans[pa]['Bmin'] = datamodel.dmarray(numpy.nan, attrs={'units':'nT'})
        ans[pa]['Bmirror'] = datamodel.dmarray(numpy.nan, attrs={'units':'nT'})

    if trace != 'LGM_CLOSED':
        return ans
        # if this is not LGM_CLOSED then don't both with any pitch angle?  true?

    #  Save field-related quantities for each Pitch Angle.
    MagEphemInfo.Pmin = Lgm_Vector.Lgm_Vector(*minB)
    MagEphemInfo.Bmin = MagEphemInfo.LstarInfo.contents.mInfo.contents.Bmin
    ans[pa]['Pmin'] = datamodel.dmarray(minB, attrs={'units':'R_E'})
    ans[pa]['Pmin'].attrs['coord_system'] = 'GSM'
    ans[pa]['Bmin'] = datamodel.dmarray(MagEphemInfo.Bmin, attrs={'units':'nT'})

    # LOOP OVER PITCH ANGLES
    for i, pa in enumerate(Alpha):

        tnow = datetime.datetime.now()
        PreStr = MagEphemInfo.LstarInfo.contents.PreStr
        PostStr = MagEphemInfo.LstarInfo.contents.PostStr
        if extended_out:
            MagEphemInfo.LstarInfo.contents.SaveShellLines = True

        # ***********************************************
        # *** not sure I fully understand this chunk, B at the mirror point*****
        # Set Pitch Angle, sin, sin^2, and Bmirror
        sa = math.sin( pa*RadPerDeg )
        sa2 = sa*sa

ans[pa]['Lstar'] = datamodel.dmarray([numpy.nan], attrs={'info':trace})
        ans[pa]['LMcIlwain'] = datamodel.dmarray([numpy.nan])
        ans[pa]['LHilton'] = datamodel.dmarray([numpy.nan])
        ans[pa]['Bmin'] = datamodel.dmarray([numpy.nan], attrs={'units':'nT'})
        ans[pa]['Bmirror'] = datamodel.dmarray([numpy.nan], attrs={'units':'nT'})
    
    if trace != 'LGM_CLOSED':
        return ans
        # if this is not LGM_CLOSED then don't both with any pitch angle?  true?

    #  Save field-related quantities for each Pitch Angle.
    MagEphemInfo.Pmin = Lgm_Vector.Lgm_Vector(*minB)
    MagEphemInfo.Bmin = MagEphemInfo.LstarInfo.contents.mInfo.contents.Bmin
    ans[pa]['Pmin'] = datamodel.dmarray(minB, attrs={'units':'R_E'})
    ans[pa]['Pmin'].attrs['coord_system'] = 'GSM'
    ans[pa]['Bmin'] = datamodel.dmarray(MagEphemInfo.Bmin, attrs={'units':'nT'})
    # LOOP OVER PITCH ANGLES
    for i, pa in enumerate(Alpha):

        tnow = datetime.datetime.now()
        PreStr = MagEphemInfo.LstarInfo.contents.PreStr
        PostStr = MagEphemInfo.LstarInfo.contents.PostStr

        # ***********************************************
        # *** not sure I fully understand this chunk, B at the mirror point*****
        # Set Pitch Angle, sin, sin^2, and Bmirror
        sa = math.sin( pa*RadPerDeg )
        sa2 = sa*sa

        # print("{0}Computing L* for Pitch Angle: Alpha[{1}] = {2} Date: {3}   UTC: {4}   Lsimple = {5:4.3}{6}\n").format(PreStr, i, MagEphemInfo.Alpha[i], date, utc, Lsimple, PostStr )
        
        if sa2!=0: #non-zero pitch angle

PreStr = MagEphemInfo.LstarInfo.contents.PreStr
        PostStr = MagEphemInfo.LstarInfo.contents.PostStr

        # ***********************************************
        # *** not sure I fully understand this chunk, B at the mirror point*****
        # Set Pitch Angle, sin, sin^2, and Bmirror
        sa = math.sin( pa*RadPerDeg )
        sa2 = sa*sa

        # print("{0}Computing L* for Pitch Angle: Alpha[{1}] = {2} Date: {3}   UTC: {4}   Lsimple = {5:4.3}{6}\n").format(PreStr, i, MagEphemInfo.Alpha[i], date, utc, Lsimple, PostStr )
        
        if sa2!=0: #non-zero pitch angle
            MagEphemInfo.LstarInfo.contents.mInfo.contents.Bm = MagEphemInfo.B/sa2
        else:
            continue
        ans[pa]['Bmirror'] = datamodel.dmarray(MagEphemInfo.LstarInfo.contents.mInfo.contents.Bm, attrs={'units':'nT'})
        ans[pa]['Bmirror'].attrs['coord_system'] = 'GSM'
        # ***********************************************
        
        # I tink this is already done
        # Lgm_Set_Coord_Transforms( Date, UTC, MagEphemInfo.LstarInfo.contents.mInfo.contents.c )

        MagEphemInfo.LstarInfo.contents.PitchAngle = pa
        MagEphemInfo.Bm[i] = MagEphemInfo.LstarInfo.contents.mInfo.contents.Bm
        # Compute L*
        if Lsimple < LstarThresh:
            Ls_vec = Lgm_Vector.Lgm_Vector(Pgsm.x, Pgsm.y, Pgsm.z)

            LS_Flag = Lgm_Lstar( pointer(Ls_vec), MagEphemInfo.LstarInfo)

            lstarinf = MagEphemInfo.LstarInfo.contents #shortcut
            #Returned from McIlwain_L but not used here

# Save nAlpha, and Alpha array to MagEphemInfo structure
    MagEphemInfo.nAlpha = len(Alpha)
    for i in range(len(Alpha)):
        MagEphemInfo.Alpha[i] = Alpha[i]

    # Set Tolerances
    Lgm_SetLstarTolerances(LstarQuality, 24, MagEphemInfo.LstarInfo )
    # *  Blocal at sat location
    MagEphemInfo.P = Pgsm

    Bvec = Lgm_Vector.Lgm_Vector()
    # Get B at the point in question
    MagEphemInfo.LstarInfo.contents.mInfo.contents.Bfield(pointer(Pgsm), pointer(Bvec),
                                                          MagEphemInfo.LstarInfo.contents.mInfo)
    ans['Bcalc'] = datamodel.dmarray(Bvec.tolist(), attrs={'units':'nT'})
    ans['Bcalc'].attrs['model'] = Bfield
    ans['Bcalc'].attrs['Kp'] = Kp
    ans['Bcalc'].attrs['coord_system'] = 'GSM'

    # save its magnitude in the structure
    MagEphemInfo.B = Bvec.magnitude()

    # check and see if the field line is closed before doing much work
    trace, northern, southern, minB, Lsimple = Closed_Field.Closed_Field(MagEphemInfo, extended_out=True)

    # presetup the ans[Angle] so that it can be filled correctly
    for pa in Alpha:
        ans[pa] = datamodel.SpaceData()
        ans[pa]['Lsimple'] = datamodel.dmarray([Lsimple])
        #sets up nans in case of Lstar failure
        ans[pa]['I'] = datamodel.dmarray([numpy.nan])