Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
def create(self, n=1):
super(Area, self).create(n)
self.oldalt[-n:] = traf.alt[-n:]
self.insdel[-n:] = False
self.insexp[-n:] = False
self.vmld = self.vm_ld*np.sqrt(self.mass/bs.traf.rho)
# summarize and convert to cas
# note: aircraft on ground may be pushed back
self.vmin = (self.phase==1)*vtas2cas(self.vmto, bs.traf.alt) + \
((self.phase==2) + (self.phase==3) + (self.phase==4))*vtas2cas(self.vmcr, bs.traf.alt) + \
(self.phase==5)*vtas2cas(self.vmld, bs.traf.alt) + (self.phase==6)*-10.0
# forwarding to tools
bs.traf.limspd, \
bs.traf.limspd_flag, \
bs.traf.limalt, \
bs.traf.limalt_flag, \
bs.traf.limvs, \
bs.traf.limvs_flag = calclimits(vtas2cas(bs.traf.pilot.tas, bs.traf.alt), \
bs.traf.gs, \
self.vmto, \
self.vmin, \
self.vmo, \
self.mmo, \
bs.traf.M, \
bs.traf.alt, \
self.hmaxact, \
bs.traf.pilot.alt, \
bs.traf.pilot.vs, \
self.maxthr, \
self.Thr_pilot, \
self.D, \
bs.traf.cas, \
self.mass, \
self.ESF, \
self.k[self.phase == ph.NA] = self.k_clean[self.phase == ph.NA]
rho = aero.vdensity(bs.traf.alt[idx_fixwing])
vtas = bs.traf.tas[idx_fixwing]
rhovs = 0.5 * rho * vtas ** 2 * self.Sref[idx_fixwing]
cl = self.mass[idx_fixwing] * aero.g0 / rhovs
self.drag[idx_fixwing] = rhovs * (
self.cd0[idx_fixwing] + self.k[idx_fixwing] * cl ** 2
)
# ----- compute maximum thrust -----
max_thrustratio_fixwing = thrust.compute_max_thr_ratio(
self.phase[idx_fixwing],
self.engbpr[idx_fixwing],
bs.traf.tas[idx_fixwing],
bs.traf.alt[idx_fixwing],
bs.traf.vs[idx_fixwing],
self.engnum[idx_fixwing] * self.engthrmax[idx_fixwing],
)
self.max_thrust[idx_fixwing] = (
max_thrustratio_fixwing
* self.engnum[idx_fixwing]
* self.engthrmax[idx_fixwing]
)
# ----- compute net thrust -----
self.thrust[idx_fixwing] = (
self.drag[idx_fixwing] + self.mass[idx_fixwing] * bs.traf.ax[idx_fixwing]
)
# ----- compute duel flow -----
thrustratio_fixwing = self.thrust[idx_fixwing] / (
# update mass
#self.mass = self.mass - self.ff*self.dt/60. # Use fuelflow in kg/min
# print bs.traf.id, self.phase, bs.traf.alt/ft, bs.traf.tas/kts, bs.traf.cas/kts, bs.traf.M, \
# self.Thr, self.D, self.ff, cl, cd, bs.traf.vs/fpm, self.ESF,self.atrans, self.maxthr, \
# self.vmto/kts, self.vmic/kts ,self.vmcr/kts, self.vmap/kts, self.vmld/kts, \
# CD0f, kf, self.hmaxact
# for aircraft on the runway and taxiways we need to know, whether they
# are prior or after their flight
self.post_flight = np.where(self.descent, True, self.post_flight)
# when landing, we would like to stop the aircraft.
bs.traf.pilot.tas = np.where((bs.traf.alt <0.5)*(self.post_flight)*self.pf_flag, 0.0, bs.traf.pilot.tas)
# the impulse for reducing the speed to 0 should only be given once,
# otherwise taxiing will be impossible afterwards
self.pf_flag = np.where ((bs.traf.alt <0.5)*(self.post_flight), False, self.pf_flag)
return
# Drag:
self.D = cd*qS
# energy share factor and crossover altitude
# conditions
epsalt = np.array([0.001]*bs.traf.ntraf)
climb = np.array(delalt > epsalt)
descent = np.array(delalt<-epsalt)
lvl = np.array(np.abs(delalt)<0.0001)*1
# energy share factor
delspd = bs.traf.pilot.tas - bs.traf.tas
selmach = bs.traf.selspd < 2.0
self.ESF = esf(bs.traf.alt, bs.traf.M, climb, descent, delspd, selmach)
# THRUST
# 1. climb: max.climb thrust in ISA conditions (p. 32, BADA User Manual 3.12)
# condition: delta altitude positive
# temperature correction for non-ISA (as soon as applied)
# ThrISA = (1-self.ctct2*(self.dtemp-self.ctct1))
# jet
# condition
cljet = np.logical_and.reduce([climb, self.jet]) * 1
# thrust
Tj = self.ctcth1* (1-(bs.traf.alt/ft)/self.ctcth2+self.ctcth3*(bs.traf.alt/ft)*(bs.traf.alt/ft))
# combine jet and default aircraft
def update(self, dt=settings.performance_dt):
"""AIRCRAFT PERFORMANCE"""
# BADA version
swbada = True
delalt = bs.traf.selalt - bs.traf.alt
# flight phase
self.phase, self.bank = phases(bs.traf.alt, bs.traf.gs, delalt,
bs.traf.cas, self.vmto, self.vmic, self.vmap, self.vmcr, self.vmld,
bs.traf.bank, bs.traf.bphase, bs.traf.swhdgsel, swbada)
# AERODYNAMICS
# Lift
qS = 0.5*bs.traf.rho*np.maximum(1.,bs.traf.tas)*np.maximum(1.,bs.traf.tas)*self.Sref
cl = self.mass*g0/(qS*np.cos(self.bank))*(self.phase!=PHASE["GD"])+ 0.*(self.phase==PHASE["GD"])
# Drag
# Drag Coefficient
# phases TO, IC, CR
cdph = self.cd0cr+self.cd2cr*(cl*cl)
# phase AP
# in case approach coefficients in OPF-Files are set to zero:
# compute drag: CD = CD0 + CDi * CL^2 and D = rho/2*VTAS^2*CD*S
self.D = cd*self.qS
# energy share factor and crossover altitude
epsalt = np.array([0.001]*bs.traf.ntraf)
self.climb = np.array(bs.traf.delalt > epsalt)
self.descent = np.array(bs.traf.delalt< -epsalt)
# crossover altitiude
bs.traf.abco = np.array(bs.traf.alt>self.atrans)
bs.traf.belco = np.array(bs.traf.alt
# Horizontal flight direction
turbhf=np.random.normal(0,self.sd[0]*timescale,bs.traf.ntraf) #[m]
# Horizontal wing direction
turbhw=np.random.normal(0,self.sd[1]*timescale,bs.traf.ntraf) #[m]
# Vertical direction
turbalt=np.random.normal(0,self.sd[2]*timescale,bs.traf.ntraf) #[m]
trkrad=np.radians(bs.traf.trk)
# Lateral, longitudinal direction
turblat=np.cos(trkrad)*turbhf-np.sin(trkrad)*turbhw #[m]
turblon=np.sin(trkrad)*turbhf+np.cos(trkrad)*turbhw #[m]
# Update the aircraft locations
bs.traf.alt = bs.traf.alt + turbalt
bs.traf.lat = bs.traf.lat + np.degrees(turblat/Rearth)
bs.traf.lon = bs.traf.lon + np.degrees(turblon/Rearth/bs.traf.coslat)
# drag coefficient
cd = self.CD0*CD0f + self.k*kf*(cl*cl)
# compute drag: CD = CD0 + CDi * CL^2 and D = rho/2*VTAS^2*CD*S
self.D = cd*self.qS
# energy share factor and crossover altitude
epsalt = np.array([0.001]*bs.traf.ntraf)
climb = np.array(delalt > epsalt)
descent = np.array(delalt< -epsalt)
# energy share factor
delspd = bs.traf.pilot.tas - bs.traf.tas
selmach = bs.traf.selspd < 2.0
self.ESF = esf(bs.traf.alt, bs.traf.M, climb, descent, delspd, selmach)
# determine thrust
self.thrust = (((bs.traf.vs*self.mass*g0)/(self.ESF*np.maximum(bs.traf.eps, bs.traf.tas))) + self.D)
# determine thrust required to fulfill requests from pilot
# self.thrust_pilot = (((bs.traf.pilot.vs*self.mass*g0)/(self.ESF*np.maximum(bs.traf.eps, bs.traf.pilot.tas))) + self.D)
self.thrust_pilot = (((bs.traf.ap.vs*self.mass*g0)/(self.ESF*np.maximum(bs.traf.eps, bs.traf.pilot.tas))) + self.D)
# maximum thrust jet (Bruenig et al., p. 66):
mt_jet = self.rated_thrust*(bs.traf.rho/rho0)**0.75
# maximum thrust prop (Raymer, p.36):
mt_prop = self.P*self.eta/np.maximum(bs.traf.eps, bs.traf.tas)
# merge
self.maxthr = mt_jet*(self.etype==1) + mt_prop*(self.etype==2)
def perf(self,simt):
if abs(simt - self.t0) >= self.dt:
self.t0 = simt
else:
return
"""Aircraft performance"""
swbada = False # no-bada version
# allocate aircraft to their flight phase
self.phase, self.bank = \
phases(bs.traf.alt, bs.traf.gs, bs.traf.delalt, \
bs.traf.cas, self.vmto, self.vmic, self.vmap, self.vmcr, self.vmld, bs.traf.bank, bs.traf.bphase, \
bs.traf.swhdgsel,swbada)
# AERODYNAMICS
# compute CL: CL = 2*m*g/(VTAS^2*rho*S)
self.qS = 0.5*bs.traf.rho*np.maximum(1.,bs.traf.tas)*np.maximum(1.,bs.traf.tas)*self.Sref
cl = self.mass*g0/(self.qS*np.cos(self.bank))*(self.phase!=6)+ 0.*(self.phase==6)
# scaling factors for CD0 and CDi during flight phases according to FAA (2005): SAGE, V. 1.5, Technical Manual
CD0f = (self.phase==1)*(self.etype==1)*coeffBS.d_CD0j[0] + \
(self.phase==2)*(self.etype==1)*coeffBS.d_CD0j[1] + \
(self.phase==3)*(self.etype==1)*coeffBS.d_CD0j[2] + \
(self.phase==4)*(self.etype==1)*coeffBS.d_CD0j[3] + \
(self.phase==5)*(self.etype==1)*(bs.traf.alt>=450)*coeffBS.d_CD0j[4] + \