diff --git a/EXOSIMS/Prototypes/Observatory.py b/EXOSIMS/Prototypes/Observatory.py index b6657594a..52c3c2abb 100644 --- a/EXOSIMS/Prototypes/Observatory.py +++ b/EXOSIMS/Prototypes/Observatory.py @@ -38,16 +38,8 @@ class Observatory(object): Path to SPK file on disk (Defaults to de432s.bsp). Attributes: - koAngleMin (astropy Quantity): - Telescope minimum keepout angle in units of deg - koAngleMinMoon (astropy Quantity): - Telescope minimum keepout angle in units of deg, for the Moon only - koAngleMinEarth (astropy Quantity): - Telescope minimum keepout angle in units of deg, for the Earth only - koAngleMax (astropy Quantity): - Telescope maximum keepout angle (for occulter) in units of deg - koAngleSmall (astropy Quantity): - Telescope keepout angle for smaller (angular size) bodies in units of deg + koAngle_SolarPanel (astropy ndarray Quantity): + Telescope minimum and maximum keepout angles (1x2 list) for solar panels in units of deg settlingTime (astropy Quantity): Instrument settling time after repoint in units of day thrust (astropy Quantity): @@ -91,9 +83,9 @@ class Observatory(object): _modtype = 'Observatory' - def __init__(self, koAngleMin=45., koAngleMinMoon=None, koAngleMinEarth=None, - koAngleMax=None, koAngleSmall=1, ko_dtStep=1, settlingTime=1, thrust=450, - slewIsp=4160., scMass=6000., dryMass=3400., coMass=5800., occulterSep=55000., skIsp=220., + def __init__(self, koAngles_SolarPanel=[0,180], + ko_dtStep=1, settlingTime=1, thrust=450, slewIsp=4160., scMass=6000., + dryMass=3400., coMass=5800., occulterSep=55000., skIsp=220., defburnPortion=0.05, constTOF=14, maxdVpct=0.02, spkpath=None, checkKeepoutEnd=True, forceStaticEphem=False, occ_dtmin=10., occ_dtmax=61., cachedir=None, **specs): @@ -108,13 +100,7 @@ def __init__(self, koAngleMin=45., koAngleMinMoon=None, koAngleMinEarth=None, assert isinstance(forceStaticEphem, bool), "forceStaticEphem must be a boolean." # default Observatory values - self.koAngleMin = float(koAngleMin)*u.deg # keepout minimum angle - koAngleMinMoon = koAngleMin if koAngleMinMoon is None else koAngleMinMoon - self.koAngleMinMoon = float(koAngleMinMoon)*u.deg # keepout minimum angle: Moon-only - koAngleMinEarth = koAngleMin if koAngleMinEarth is None else koAngleMinEarth - self.koAngleMinEarth = float(koAngleMinEarth)*u.deg# keepout minimum angle: Earth-only - self.koAngleMax = float(koAngleMax)*u.deg if koAngleMax is not None else koAngleMax # keepout maximum angle (occulter) - self.koAngleSmall = float(koAngleSmall)*u.deg # keepout angle for smaller bodies + self.koAngles_SolarPanel = [float(x) for x in koAngles_SolarPanel]*u.deg #solar panel keepout angles self.ko_dtStep = float(ko_dtStep)*u.d # time step for generating koMap of stars (day) self.settlingTime = float(settlingTime)*u.d # instru. settling time after repoint self.thrust = float(thrust)*u.mN # occulter slew thrust (mN) @@ -124,10 +110,10 @@ def __init__(self, koAngleMin=45., koAngleMinMoon=None, koAngleMinEarth=None, self.coMass = float(coMass)*u.kg # telescope mass (kg) self.occulterSep = float(occulterSep)*u.km # occulter-telescope distance (km) self.skIsp = float(skIsp)*u.s # station-keeping Isp (s) - self.defburnPortion = float(defburnPortion) # default burn portion - self.checkKeepoutEnd = bool(checkKeepoutEnd)# true if keepout called at obs end - self.forceStaticEphem = bool(forceStaticEphem)# boolean used to force static ephem - self.constTOF = np.array(constTOF,ndmin=1)*u.d # starshade constant slew time (days) + self.defburnPortion = float(defburnPortion) # default burn portion + self.checkKeepoutEnd = bool(checkKeepoutEnd) # true if keepout called at obs end + self.forceStaticEphem = bool(forceStaticEphem) # boolean used to force static ephem + self.constTOF = np.array(constTOF,ndmin=1)*u.d # starshade constant slew time (days) self.occ_dtmin = float(occ_dtmin)*u.d # Minimum occulter slew time (days) self.occ_dtmax = float(occ_dtmax)*u.d # Maximum occulter slew time (days) self.maxdVpct = float(maxdVpct) # Maximum deltaV percent @@ -404,7 +390,7 @@ def orbit(self, currentTime, eclip=False): return r_obs - def keepout(self, TL, sInds, currentTime, returnExtra=False): + def keepout(self, TL, sInds, currentTime, koangles, returnExtra=False): """Finds keepout Boolean values for stars of interest. This method returns the keepout Boolean values for stars of interest, where @@ -417,27 +403,32 @@ def keepout(self, TL, sInds, currentTime, returnExtra=False): Integer indices of the stars of interest currentTime (astropy Time array): Current absolute mission time in MJD + koangles (astropy Quantity ndarray): + s x 4 x 2 array where s is the number of starlight suppression systems as + defined in the Optical System. Each of the remaining 4 x 2 arrays are system + specific koAngles for the Sun, Earth, Moon, and small bodies (4), each with a + minimum and maximum value (2) in units of deg. returnExtra (boolean): Optional flag, default False, set True to return additional information: r_body (astropy Quantity array): - 11 x n x 3 array where n is len(currentTime) of heliocentric + 11 x m x 3 array where m is len(currentTime) of heliocentric equatorial Cartesian elements of the Sun, Moon, Earth and Mercury->Pluto r_targ (astropy Quantity array): - n x 3 array where n is len(currentTime) or 1 if staticStars is true in - TargetList of heliocentric equatorial Cartesian coords of target + m x n x 3 array where m is len(currentTime) or 1 if staticStars is true in + TargetList of heliocentric equatorial Cartesian coords of target and n is the + len(sInds) culprit (float ndarray): - m x n x 11 array of boolean integer values identifying which body + s x n x m x 12 array of boolean integer values identifying which body is responsible for keepout (when equal to 1). m is number of targets - and n is len(currentTime). Last dimension is ordered same as r_body - koangles (astropy quantity ndarray): - 11 element array of keepouts used for each body. Same ordering as - r_body. + and n is len(currentTime). Last dimension is ordered same as r_body, with + an extra line for solar panels being the culprit + koangleArray (astropy quantity ndarray): + s x 11 x 2 element array of minimum and maximum keepouts used for each body. + Same ordering as r_body. Returns: boolean ndarray: - True is a target unobstructed and observable, and False is a - target unobservable due to obstructions in the keepout zone. - - Note: If multiple times and targets, currentTime and sInds sizes must match. + s x n x m array of boolean values. True is a target unobstructed and observable, + and False is a target unobservable due to obstructions in the keepout zone. """ @@ -445,77 +436,80 @@ def keepout(self, TL, sInds, currentTime, returnExtra=False): if currentTime.size > 1: if np.all(currentTime == currentTime[0]): currentTime = currentTime[0] - + # cast sInds to array sInds = np.array(sInds, ndmin=1, copy=False) # get all array sizes nStars = sInds.size nTimes = currentTime.size + nSystems = koangles.shape[0] nBodies = 11 - assert nStars == 1 or nTimes == 1 or nTimes == nStars, \ - "If multiple times and targets, currentTime and sInds sizes must match" - + # observatory positions vector in heliocentric equatorial frame - r_obs = self.orbit(currentTime) + r_obs = self.orbit(currentTime) # (m x 3) # traget star positions vector in heliocentric equatorial frame - r_targ = TL.starprop(sInds, currentTime) + r_targ = TL.starprop(sInds, currentTime).reshape(nTimes,nStars,3) # (m x n x 3) # body positions vector in heliocentric equatorial frame r_body = np.array([ - self.solarSystem_body_position(currentTime, 'Sun').to('AU').value, - self.solarSystem_body_position(currentTime, 'Moon').to('AU').value, - self.solarSystem_body_position(currentTime, 'Earth').to('AU').value, - self.solarSystem_body_position(currentTime, 'Mercury').to('AU').value, - self.solarSystem_body_position(currentTime, 'Venus').to('AU').value, - self.solarSystem_body_position(currentTime, 'Mars').to('AU').value, - self.solarSystem_body_position(currentTime, 'Jupiter').to('AU').value, - self.solarSystem_body_position(currentTime, 'Saturn').to('AU').value, - self.solarSystem_body_position(currentTime, 'Uranus').to('AU').value, - self.solarSystem_body_position(currentTime, 'Neptune').to('AU').value, - self.solarSystem_body_position(currentTime, 'Pluto').to('AU').value])*u.AU + self.solarSystem_body_position(currentTime, 'Sun').to('AU').value, + self.solarSystem_body_position(currentTime, 'Moon').to('AU').value, + self.solarSystem_body_position(currentTime, 'Earth').to('AU').value, + self.solarSystem_body_position(currentTime, 'Mercury').to('AU').value, + self.solarSystem_body_position(currentTime, 'Venus').to('AU').value, + self.solarSystem_body_position(currentTime, 'Mars').to('AU').value, + self.solarSystem_body_position(currentTime, 'Jupiter').to('AU').value, + self.solarSystem_body_position(currentTime, 'Saturn').to('AU').value, + self.solarSystem_body_position(currentTime, 'Uranus').to('AU').value, + self.solarSystem_body_position(currentTime, 'Neptune').to('AU').value, + self.solarSystem_body_position(currentTime, 'Pluto').to('AU').value])*u.AU # position vectors wrt spacecraft - r_targ = (r_targ - r_obs).to('pc') - r_body = (r_body - r_obs).to('AU') + r_targ = (r_targ - r_obs.reshape(nTimes, 1,3) ).to('pc') # (m x n x 3) + r_body = (r_body - r_obs.reshape( 1,nTimes,3) ).to('AU') # (11 x m x 3) # unit vectors wrt spacecraft - u_targ = (r_targ.value.T/np.linalg.norm(r_targ, axis=-1)).T - u_body = (r_body.value.T/np.linalg.norm(r_body, axis=-1).T).T - - # create koangles for all bodies, set by telescope minimum keepout angle for - # brighter objects (Sun, Moon, Earth) and defaults to 1 degree for other bodies - koangles = np.ones(nBodies)*self.koAngleMin - # allow Moon, Earth to be set individually (default to koAngleMin) - koangles[1] = self.koAngleMinMoon - koangles[2] = self.koAngleMinEarth - # keepout angle for small bodies (other planets) - koangles[3:] = self.koAngleSmall + u_targ = (r_targ.value/np.linalg.norm(r_targ, axis=-1, keepdims=True)) + u_body = (r_body.value/np.linalg.norm(r_body, axis=-1, keepdims=True)) + + # create array of koangles for all bodies, using minimum and maximum keepout + # angles of each starlight suppression system in the telescope for + # bright objects (Sun, Moon, Earth, other small bodies) + koangleArray = np.zeros([nSystems, nBodies, 2]) + koangleArray[:,0:3 ,:] = koangles[:,0:3,:] + koangleArray[:,3: ,:] = koangles[:,3,:].reshape(nSystems,1,2) #small bodies have same values + koangleArray = koangleArray*u.deg # find angles and make angle comparisons to build kogood array: # if bright objects have an angle with the target vector less than koangle # (e.g. pi/4) they are in the field of view and the target star may not be - # observed, thus ko associated with this target becomes False - nkogood = np.maximum(nStars, nTimes) - kogood = np.array([True]*nkogood) - culprit = np.zeros([nkogood, nBodies]) - for i in xrange(nkogood): - u_b = u_body[:,0,:] if nTimes == 1 else u_body[:,i,:] - u_t = u_targ[0,:] if nStars == 1 else u_targ[i,:] - angles = np.arccos(np.clip(np.dot(u_b, u_t), -1, 1))*u.rad - culprit[i,:] = (angles < koangles) - # if this mode has an occulter, check maximum keepout angle for the Sun - if self.koAngleMax is not None: - culprit[i,0] = (culprit[i,0] or (angles[0] > self.koAngleMax)) - if np.any(culprit[i,:]): - kogood[i] = False - - # check to make sure all elements in kogood are Boolean - trues = [isinstance(element, np.bool_) for element in kogood] + # observed, thus ko associated with this target becomes False. + kogood = np.ones( [nSystems, nStars, nTimes], dtype=bool) # (s x n x m) + culprit = np.zeros([nSystems, nStars, nTimes, nBodies+1]) # (s x n x m x 12) + # running loop for nSystems, nStars, and nTimes (three loops total) + for s in np.arange(nSystems): + for n in np.arange(nStars): + for m in np.arange(nTimes): + # unit vectors for the 11 bodies and the nth target at the mth time + u_b = u_body[:,m,:] + u_t = u_targ[m,n,:] + # relative angle between the target and bright body look vectors + angles = np.arccos(np.clip(np.dot(u_b, u_t), -1, 1))*u.rad + # create array of "culprits" that prevent a target from being observed + culprit[s,n,m,:-1] = (angleskoangleArray[s,:,1]) + # adding solar panel restrictions as a final culprit + culprit[s,n,m,-1] = (angles[0]self.koAngles_SolarPanel[1]) + # if any bright body obstructs, kogood becomes False + if np.any(culprit[s,n,m,:]): + kogood[s,n,m] = False + + # checking that all ko elements are boolean + trues = [isinstance(element, np.bool_) for element in kogood.flatten()] assert all(trues), "An element of kogood is not Boolean" - + if returnExtra: - return kogood, r_body, r_targ, culprit, koangles + return kogood, r_body, r_targ, culprit, koangleArray else: return kogood - def generate_koMap(self,TL,missionStart,missionFinishAbs): + def generate_koMap(self,TL,missionStart,missionFinishAbs,koangles): """Creates keepout map for all targets throughout mission lifetime. This method returns a binary map showing when all stars in the given @@ -529,8 +523,11 @@ def generate_koMap(self,TL,missionStart,missionFinishAbs): Absolute start of mission time in MJD missionFinishAbs (astropy Time array): Absolute end of mission time in MJD - mode (dict): - Selected observing mode + koangles (astropy Quantity ndarray): + s x 4 x 2 array where s is the number of starlight suppression systems as + defined in the Optical System. Each of the remaining 4 x 2 arrays are system + specific koAngles for the Sun, Earth, Moon, and small bodies (4), each with a + minimum and maximum value (2) in units of deg. Returns: tuple: @@ -543,13 +540,14 @@ def generate_koMap(self,TL,missionStart,missionFinishAbs): """ # generating hash name filename = 'koMap_' - atts = ['koAngleMin', 'koAngleMinMoon', 'koAngleMinEarth', 'koAngleMax', 'koAngleSmall', 'ko_dtStep'] + atts = ['koAngles_SolarPanel','ko_dtStep'] extstr = '' for att in sorted(atts, key=str.lower): if not callable(getattr(self, att)): extstr += '%s: ' % att + str(getattr(self, att)) + ' ' extstr += '%s: ' % 'missionStart' + str(missionStart) + ' ' extstr += '%s: ' % 'missionFinishAbs' + str(missionFinishAbs) + ' ' + extstr += '%s: ' % 'koangles' + str(koangles) + ' ' extstr += '%s: ' % 'Name' + str(getattr(TL, 'Name')) + ' ' ext = hashlib.md5(extstr.encode('utf-8')).hexdigest() filename += ext @@ -574,10 +572,7 @@ def generate_koMap(self,TL,missionStart,missionFinishAbs): self.vprint('Cached keepout map file not found at "%s".' % koPath) # looping over all stars to generate map of when all stars are observable self.vprint('Starting keepout calculations for %s stars.' % TL.nStars) - koMap = np.zeros([TL.nStars,len(koTimes)]) - for n in range(TL.nStars): - koMap[n,:] = self.keepout(TL,n,koTimes,False) - if not n % 50: self.vprint(' [%s / %s] completed.' % (n,TL.nStars)) + koMap = self.keepout(TL,np.arange(TL.nStars),koTimes,koangles,False) A = {'koMap':koMap} with open(koPath, 'wb') as f: pickle.dump(A, f) @@ -586,7 +581,7 @@ def generate_koMap(self,TL,missionStart,missionFinishAbs): return koMap,koTimes - def calculate_observableTimes(self, TL, sInds, currentTime, koMap, koTimes, mode): + def calculate_observableTimes(self, TL, sInds, currentTime, koMaps, koTimes, mode): """Returns the next window of time during which targets are observable This method returns a nx2 ndarray of times for every star given in the @@ -601,8 +596,11 @@ def calculate_observableTimes(self, TL, sInds, currentTime, koMap, koTimes, mode Integer indices of the stars of interest currentTime (astropy Time array): Current absolute mission time in MJD - koMap (integer ndarray nxm): - Keepout values for n stars throughout time range of length m + koMaps (dict): + Keepout values for n stars throughout time range of length m, + key names being the system names specified in mode. + True is a target unobstructed and observable, and False is a + target unobservable due to obstructions in the keepout zone. koTimes (astropy Time ndarray): Absolute MJD mission times from start to end in steps of 1 d mode (dict): @@ -619,6 +617,10 @@ def calculate_observableTimes(self, TL, sInds, currentTime, koMap, koTimes, mode else: nextObTimes = np.ones(len(sInds))*currentTime.value nextObTimes = Time(nextObTimes,format='mjd',scale='tai') #converting to astropy MJD time array + #find appropriate koMap + systName = mode['syst']['name'] + koMap = koMaps[systName] + # finding observable times observableTimes = self.find_nextObsWindow(TL,sInds,nextObTimes,koMap,koTimes).value observableTimesNorm = observableTimes - nextObTimes.value #days since currentTime diff --git a/EXOSIMS/Prototypes/OpticalSystem.py b/EXOSIMS/Prototypes/OpticalSystem.py index 17e794beb..3b180eaa2 100644 --- a/EXOSIMS/Prototypes/OpticalSystem.py +++ b/EXOSIMS/Prototypes/OpticalSystem.py @@ -60,6 +60,14 @@ class OpticalSystem(object): Mission observing modes attributes cachedir (str): Path to EXOSIMS cache directory + koAngles_Sun (astropy Quantity): + Telescope minimum and maximum keepout angle in units of deg + koAngles_Earth (astropy Quantity): + Telescope minimum and maximum keepout angle in units of deg, for the Earth only + koAngles_Moon (astropy Quantity): + Telescope minimum and maximum keepout angle in units of deg, for the Moon only + koAngles_Small (astropy Quantity): + Telescope minimum and maximum keepout angle (for small bodies) in units of deg Common science instrument attributes: name (string): @@ -196,6 +204,7 @@ def __init__(self, obscurFac=0.1, shapeFac=np.pi/4, pupilDiam=4, intCutoff=50, core_thruput=0.1, core_contrast=1e-10, core_platescale=None, PSF=np.ones((3,3)), ohTime=1, observingModes=None, SNR=5, timeMultiplier=1., IWA=None, OWA=None, ref_dMag=3, ref_Time=0, cachedir=None, + koAngles_Sun=[0,180], koAngles_Earth=[0,180], koAngles_Moon=[0,180], koAngles_Small=[0,180], use_char_minintTime=False, **specs): #start the outspec @@ -334,6 +343,12 @@ def __init__(self, obscurFac=0.1, shapeFac=np.pi/4, pupilDiam=4, intCutoff=50, if nsyst == 0: lam, BW = syst.get('lam').value, syst.get('BW') + # get keepout angles for specific instrument + syst['koAngles_Sun'] = [float(x) for x in syst.get('koAngles_Sun', koAngles_Sun)]*u.deg + syst['koAngles_Earth'] = [float(x) for x in syst.get('koAngles_Earth',koAngles_Earth)]*u.deg + syst['koAngles_Moon'] = [float(x) for x in syst.get('koAngles_Moon', koAngles_Moon)]*u.deg + syst['koAngles_Small'] = [float(x) for x in syst.get('koAngles_Small',koAngles_Small)]*u.deg + # get coronagraph input parameters syst = self.get_coro_param(syst, 'occ_trans') syst = self.get_coro_param(syst, 'core_thruput') diff --git a/EXOSIMS/Prototypes/SurveySimulation.py b/EXOSIMS/Prototypes/SurveySimulation.py index 018f01ba6..ca0307351 100644 --- a/EXOSIMS/Prototypes/SurveySimulation.py +++ b/EXOSIMS/Prototypes/SurveySimulation.py @@ -290,16 +290,34 @@ def __init__(self, scriptfile=None, ntFlux=1, nVisitsMax=5, charMargin=0.15, #Generate File Hashnames and loction self.cachefname = self.generateHashfName(specs) - # getting keepout map for entire mission - startTime = self.TimeKeeping.missionStart.copy() - endTime = self.TimeKeeping.missionFinishAbs.copy() - if not(nokoMap): - self.koMap,self.koTimes = self.Observatory.generate_koMap(TL,startTime,endTime) - # choose observing modes selected for detection (default marked with a flag) allModes = OS.observingModes det_mode = list(filter(lambda mode: mode['detectionMode'] == True, allModes))[0] self.mode = det_mode + + # getting keepout map for entire mission + startTime = self.TimeKeeping.missionStart.copy() + endTime = self.TimeKeeping.missionFinishAbs.copy() + + nSystems = len(allModes) + systNames = np.unique([allModes[x]['syst']['name'] for x in np.arange(nSystems)]).tolist() + koStr = list(filter(lambda syst: syst.startswith('koAngles_') , allModes[0]['syst'].keys())) + koangles = np.zeros([len(systNames),4,2]) + tmpNames = list(systNames) + cnt = 0 + + for x in np.arange(nSystems): + name = allModes[x]['syst']['name'] + if name in tmpNames: + koangles[cnt] = np.asarray([allModes[x]['syst'][k] for k in koStr]) + cnt += 1 + tmpNames.remove(name) + + if not(nokoMap): + koMaps,self.koTimes = self.Observatory.generate_koMap(TL,startTime,endTime,koangles) + self.koMaps = {} + for x,n in enumerate(systNames): + self.koMaps[n] = koMaps[x,:,:] # Precalculating intTimeFilter sInds = np.arange(TL.nStars) #Initialize some sInds array @@ -471,7 +489,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -550,7 +568,9 @@ def next_target(self, old_sInd, mode): # allocate settling time + overhead time tmpCurrentTimeAbs = TK.currentTimeAbs.copy() + Obs.settlingTime + mode['syst']['ohTime'] tmpCurrentTimeNorm = TK.currentTimeNorm.copy() + Obs.settlingTime + mode['syst']['ohTime'] - + + #create appropriate koMap + koMap = self.koMaps[mode['syst']['name']] # look for available targets # 1. initialize arrays @@ -566,7 +586,7 @@ def next_target(self, old_sInd, mode): sd = None if OS.haveOcculter == True: sd = Obs.star_angularSep(TL, old_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMap,self.koTimes,mode) + obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMaps,self.koTimes,mode) slewTimes = Obs.calculate_slewTimes(TL, old_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -582,7 +602,7 @@ def next_target(self, old_sInd, mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(startTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool except:#If there are no target stars to observe @@ -603,7 +623,7 @@ def next_target(self, old_sInd, mode): else: intTimes[sInds] = self.calc_targ_intTime(sInds, startTimes[sInds], mode) sInds = sInds[np.where(intTimes[sInds] <= maxIntTime)] # Filters targets exceeding end of OB - endTimes = startTimes + intTimes + endTimes = tmpCurrentTimeAbs.copy() + intTimes if maxIntTime.value <= 0: sInds = np.asarray([],dtype=int) @@ -617,7 +637,7 @@ def next_target(self, old_sInd, mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(endTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool except: @@ -1420,8 +1440,11 @@ def observation_characterization(self, sInd, mode): startTime = TK.currentTimeAbs.copy() + mode['syst']['ohTime'] + Obs.settlingTime startTimeNorm = TK.currentTimeNorm.copy() + mode['syst']['ohTime'] + Obs.settlingTime # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime) - + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + koMap = self.koMaps[mode['syst']['name']] + tochar[tochar] = koMap[sInd][koTimeInd] + # 2/ if any planet to characterize, find the characterization times # at the detected fEZ, dMag, and WA if np.any(tochar): @@ -1441,9 +1464,16 @@ def observation_characterization(self, sInd, mode): # planets to characterize tochar = ((totTimes > 0) & (totTimes <= OS.intCutoff) & (endTimesNorm <= TK.OBendTimes[TK.OBnumber])) + # 3/ is target still observable at the end of any char time? if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar]) + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] # 4/ if yes, allocate the overhead time, and perform the characterization # for the maximum char time diff --git a/EXOSIMS/Prototypes/TargetList.py b/EXOSIMS/Prototypes/TargetList.py index b55f15918..ef6130fb6 100644 --- a/EXOSIMS/Prototypes/TargetList.py +++ b/EXOSIMS/Prototypes/TargetList.py @@ -666,27 +666,36 @@ def starprop(self, sInds, currentTime, eclip=False): False, corresponding to heliocentric equatorial frame. Returns: - astropy Quantity nx3 array: + r_targ (astropy Quantity array): Target star positions vector in heliocentric equatorial (default) - or ecliptic frame in units of pc + or ecliptic frame in units of pc. Will return an m x n x 3 array + where m is size of currentTime, n is size of sInds. If either m or + n is 1, will return n x 3 or m x 3. Note: Use eclip=True to get ecliptic coordinates. """ + # if multiple time values, check they are different otherwise reduce to scalar + if currentTime.size > 1: + if np.all(currentTime == currentTime[0]): + currentTime = currentTime[0] + # cast sInds to array sInds = np.array(sInds, ndmin=1, copy=False) - # if the starprop_static method was created (staticStars is True), then use it - if self.starprop_static is not None: - return self.starprop_static(sInds, currentTime, eclip) - # get all array sizes nStars = sInds.size nTimes = currentTime.size - assert nStars==1 or nTimes==1 or nTimes==nStars, \ - "If multiple times and targets, currentTime and sInds sizes must match" - + + # if the starprop_static method was created (staticStars is True), then use it + if self.starprop_static is not None: + r_targ = self.starprop_static(sInds, currentTime, eclip) + if (nTimes == 1 or nStars == 1 or nTimes == nStars): + return r_targ + else: + return np.tile(r_targ, (nTimes, 1, 1)) + # target star ICRS coordinates coord_old = self.coords[sInds] # right ascension and declination @@ -702,17 +711,34 @@ def starprop(self, sInds, currentTime, eclip=False): v = mu0/self.parx[sInds]*u.AU + r0*self.rv[sInds] # set J2000 epoch j2000 = Time(2000., format='jyear') - # target star positions vector in heliocentric equatorial frame - dr = v*(currentTime.mjd - j2000.mjd)*u.day - r_targ = (coord_old.cartesian.xyz + dr).T.to('pc') - - if eclip: - # transform to heliocentric true ecliptic frame - coord_new = SkyCoord(r_targ[:,0], r_targ[:,1], r_targ[:,2], - representation='cartesian') - r_targ = coord_new.heliocentrictrueecliptic.cartesian.xyz.T.to('pc') - - return r_targ + + # if only 1 time in currentTime + if (nTimes == 1 or nStars == 1 or nTimes == nStars): + # target star positions vector in heliocentric equatorial frame + dr = v*(currentTime.mjd - j2000.mjd)*u.day + r_targ = (coord_old.cartesian.xyz + dr).T.to('pc') + + if eclip: + # transform to heliocentric true ecliptic frame + coord_new = SkyCoord(r_targ[:,0], r_targ[:,1], r_targ[:,2], + representation='cartesian') + r_targ = coord_new.heliocentrictrueecliptic.cartesian.xyz.T.to('pc') + return r_targ + + # create multi-dimensional array for r_targ + else: + # target star positions vector in heliocentric equatorial frame + r_targ = np.zeros([nTimes,nStars,3])*u.pc + for i,m in enumerate(currentTime): + dr = v*(m.mjd - j2000.mjd)*u.day + r_targ[i,:,:] = (coord_old.cartesian.xyz + dr).T.to('pc') + + if eclip: + # transform to heliocentric true ecliptic frame + coord_new = SkyCoord(r_targ[i,:,0], r_targ[i,:,1], r_targ[i,:,2], + representation='cartesian') + r_targ[i,:,:] = coord_new.heliocentrictrueecliptic.cartesian.xyz.T.to('pc') + return r_targ def starMag(self, sInds, lam): """Calculates star visual magnitudes with B-V color using empirical fit diff --git a/EXOSIMS/SurveySimulation/ExoC_Scheduler.py b/EXOSIMS/SurveySimulation/ExoC_Scheduler.py index 834bf715e..31e2e01ab 100644 --- a/EXOSIMS/SurveySimulation/ExoC_Scheduler.py +++ b/EXOSIMS/SurveySimulation/ExoC_Scheduler.py @@ -173,7 +173,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -271,7 +271,10 @@ def observation_characterization(self, sInd, mode, mode_index): startTime = TK.currentTimeAbs.copy() + mode['syst']['ohTime'] + Obs.settlingTime startTimeNorm = TK.currentTimeNorm.copy() + mode['syst']['ohTime'] + Obs.settlingTime # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + koMap = self.koMaps[mode['syst']['name']] + tochar[tochar] = koMap[sInd][koTimeInd] # 2/ if any planet to characterize, find the characterization times # at the detected fEZ, dMag, and WA @@ -294,7 +297,13 @@ def observation_characterization(self, sInd, mode, mode_index): (endTimesNorm <= TK.OBendTimes[TK.OBnumber])) # 3/ is target still observable at the end of any char time? if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar]) + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] # 4/ if yes, allocate the overhead time, and perform the characterization # for the maximum char time diff --git a/EXOSIMS/SurveySimulation/SS_char_only.py b/EXOSIMS/SurveySimulation/SS_char_only.py index 585c97118..dd5dc9a5b 100644 --- a/EXOSIMS/SurveySimulation/SS_char_only.py +++ b/EXOSIMS/SurveySimulation/SS_char_only.py @@ -286,6 +286,9 @@ def observation_characterization(self, sInd, mode): Obs = self.Observatory TK = self.TimeKeeping + # selecting appropriate koMap + koMap = self.koMaps[mode['syst']['name']] + # find indices of planets around the target pInds = np.where(SU.plan2star == sInd)[0] # get the last detected planets, and check if there was a FA @@ -332,7 +335,10 @@ def observation_characterization(self, sInd, mode): startTime = TK.currentTimeAbs startTimeNorm = TK.currentTimeNorm # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime, mode) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + tochar[tochar] = koMap[sInd][koTimeInd] + # 2/ if any planet to characterize, find the characterization times if np.any(tochar): @@ -361,7 +367,13 @@ def observation_characterization(self, sInd, mode): # 3/ is target still observable at the end of any char time? if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar], mode) + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] # 4/ if yes, perform the characterization for the maximum char time if np.any(tochar): diff --git a/EXOSIMS/SurveySimulation/SS_char_only2.py b/EXOSIMS/SurveySimulation/SS_char_only2.py index 08b42b3e8..a2462c43d 100644 --- a/EXOSIMS/SurveySimulation/SS_char_only2.py +++ b/EXOSIMS/SurveySimulation/SS_char_only2.py @@ -227,6 +227,9 @@ def observation_characterization(self, sInd, mode): Obs = self.Observatory TK = self.TimeKeeping + # selecting appropriate koMap + koMap = self.koMaps[mode['syst']['name']] + # find indices of planets around the target pInds = np.where(SU.plan2star == sInd)[0] # get the last detected planets, and check if there was a FA @@ -273,7 +276,9 @@ def observation_characterization(self, sInd, mode): startTime = TK.currentTimeAbs startTimeNorm = TK.currentTimeNorm # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime, mode) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + tochar[tochar] = koMap[sInd][koTimeInd] # 2/ if any planet to characterize, find the characterization times if np.any(tochar): @@ -302,7 +307,13 @@ def observation_characterization(self, sInd, mode): # 3/ is target still observable at the end of any char time? if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar], mode) + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] # 4/ if yes, perform the characterization for the maximum char time if np.any(tochar): diff --git a/EXOSIMS/SurveySimulation/linearJScheduler.py b/EXOSIMS/SurveySimulation/linearJScheduler.py index e6f13fded..2c864dbde 100644 --- a/EXOSIMS/SurveySimulation/linearJScheduler.py +++ b/EXOSIMS/SurveySimulation/linearJScheduler.py @@ -83,6 +83,9 @@ def next_target(self, old_sInd, mode): # create DRM DRM = {} + #create appropriate koMap + koMap = self.koMaps[mode['syst']['name']] + # allocate settling time + overhead time tmpCurrentTimeAbs = TK.currentTimeAbs.copy() + Obs.settlingTime + mode['syst']['ohTime'] tmpCurrentTimeNorm = TK.currentTimeNorm.copy() + Obs.settlingTime + mode['syst']['ohTime'] @@ -101,7 +104,7 @@ def next_target(self, old_sInd, mode): sd = None if OS.haveOcculter == True: sd = Obs.star_angularSep(TL, old_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMap,self.koTimes,mode) + obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMaps,self.koTimes,mode) slewTimes = Obs.calculate_slewTimes(TL, old_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -117,7 +120,7 @@ def next_target(self, old_sInd, mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(startTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool #DELETE koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] @@ -140,7 +143,7 @@ def next_target(self, old_sInd, mode): # else: intTimes[sInds] = self.calc_targ_intTime(sInds, startTimes[sInds], mode) sInds = sInds[np.where(intTimes[sInds] <= maxIntTime)] # Filters targets exceeding end of OB - endTimes = startTimes + intTimes + endTimes = tmpCurrentTimeAbs.copy() + intTimes if maxIntTime.value <= 0: sInds = np.asarray([],dtype=int) @@ -154,7 +157,7 @@ def next_target(self, old_sInd, mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(endTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool #DELETE koTimeInd = np.where(np.round(endTimes[0].value)-self.koTimes.value==0)[0][0]#koTimeInd[0][0] # find indice where koTime is endTime[0] @@ -460,7 +463,10 @@ def observation_characterization(self, sInd, mode): startTime = TK.currentTimeAbs.copy() + mode['syst']['ohTime'] + Obs.settlingTime startTimeNorm = TK.currentTimeNorm.copy() + mode['syst']['ohTime'] + Obs.settlingTime # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + koMap = self.koMaps[mode['syst']['name']] + tochar[tochar] = koMap[sInd][koTimeInd] # 2/ if any planet to characterize, find the characterization times # at the detected fEZ, dMag, and WA @@ -488,7 +494,13 @@ def observation_characterization(self, sInd, mode): (endTimesNorm <= TK.OBendTimes[TK.OBnumber])) # 3/ is target still observable at the end of any char time? if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar]) + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] # 4/ if yes, allocate the overhead time, and perform the characterization # for the maximum char time diff --git a/EXOSIMS/SurveySimulation/linearJScheduler_3DDPC.py b/EXOSIMS/SurveySimulation/linearJScheduler_3DDPC.py index 6ca7d1408..309d7d76d 100644 --- a/EXOSIMS/SurveySimulation/linearJScheduler_3DDPC.py +++ b/EXOSIMS/SurveySimulation/linearJScheduler_3DDPC.py @@ -164,7 +164,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL, np.arange(TL.nStars), startTimes, self.koMap, self.koTimes, self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL, np.arange(TL.nStars), startTimes, self.koMaps, self.koTimes, self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -234,6 +234,9 @@ def next_target(self, old_sInd, modes): # create DRM DRM = {} + # selecting appropriate koMap + koMap = self.koMaps[modes[0]['syst']['name']] + # allocate settling time + overhead time tmpCurrentTimeAbs = TK.currentTimeAbs.copy() + Obs.settlingTime + modes[0]['syst']['ohTime'] tmpCurrentTimeNorm = TK.currentTimeNorm.copy() + Obs.settlingTime + modes[0]['syst']['ohTime'] @@ -256,7 +259,7 @@ def next_target(self, old_sInd, modes): sd = None if OS.haveOcculter == True: sd = Obs.star_angularSep(TL, old_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMap, self.koTimes, mode) + obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMaps, self.koTimes, mode) slewTimes = Obs.calculate_slewTimes(TL, old_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -270,7 +273,7 @@ def next_target(self, old_sInd, modes): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - mode_sInds = mode_sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[mode_sInds])[0]]# filters inds by koMap #verified against v1.35 + mode_sInds = mode_sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[mode_sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe mode_sInds = np.asarray([],dtype=int) @@ -305,7 +308,7 @@ def next_target(self, old_sInd, modes): if len(mode_sInds.tolist()) > 0 and Obs.checkKeepoutEnd: try: # endTimes may exist past koTimes so we have an exception to hand this case koTimeInd = np.where(np.round(endTimes[0].value)-self.koTimes.value==0)[0][0]#koTimeInd[0][0] # find indice where koTime is endTime[0] - mode_sInds = mode_sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[mode_sInds])[0]]# filters inds by koMap #verified against v1.35 + mode_sInds = mode_sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[mode_sInds])[0]]# filters inds by koMap #verified against v1.35 except: mode_sInds = np.asarray([],dtype=int) diff --git a/EXOSIMS/SurveySimulation/linearJScheduler_3DDPC_sotoSS.py b/EXOSIMS/SurveySimulation/linearJScheduler_3DDPC_sotoSS.py index 421501a9d..57f43550d 100644 --- a/EXOSIMS/SurveySimulation/linearJScheduler_3DDPC_sotoSS.py +++ b/EXOSIMS/SurveySimulation/linearJScheduler_3DDPC_sotoSS.py @@ -166,7 +166,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL, np.arange(TL.nStars), startTimes, self.koMap, self.koTimes, self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL, np.arange(TL.nStars), startTimes, self.koMaps, self.koTimes, self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -236,6 +236,9 @@ def next_target(self, old_sInd, modes): # create DRM DRM = {} + # selecting appropriate koMap + koMap = self.koMaps[char_mode['syst']['name']] + # allocate settling time + overhead time tmpCurrentTimeAbs = TK.currentTimeAbs.copy() + Obs.settlingTime + modes[0]['syst']['ohTime'] tmpCurrentTimeNorm = TK.currentTimeNorm.copy() + Obs.settlingTime + modes[0]['syst']['ohTime'] @@ -258,7 +261,7 @@ def next_target(self, old_sInd, modes): sd = None if OS.haveOcculter == True: sd = Obs.star_angularSep(TL, old_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMap, self.koTimes, mode) + obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMaps, self.koTimes, mode) slewTimes = Obs.calculate_slewTimes(TL, old_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -272,7 +275,7 @@ def next_target(self, old_sInd, modes): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - mode_sInds = mode_sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[mode_sInds])[0]]# filters inds by koMap #verified against v1.35 + mode_sInds = mode_sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[mode_sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe mode_sInds = np.asarray([],dtype=int) @@ -307,7 +310,7 @@ def next_target(self, old_sInd, modes): if len(mode_sInds.tolist()) > 0 and Obs.checkKeepoutEnd: try: # endTimes may exist past koTimes so we have an exception to hand this case koTimeInd = np.where(np.round(endTimes[0].value)-self.koTimes.value==0)[0][0]#koTimeInd[0][0] # find indice where koTime is endTime[0] - mode_sInds = mode_sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[mode_sInds])[0]]# filters inds by koMap #verified against v1.35 + mode_sInds = mode_sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[mode_sInds])[0]]# filters inds by koMap #verified against v1.35 except: mode_sInds = np.asarray([],dtype=int) diff --git a/EXOSIMS/SurveySimulation/linearJScheduler_DDPC.py b/EXOSIMS/SurveySimulation/linearJScheduler_DDPC.py index 54aa10817..ae7c83fec 100644 --- a/EXOSIMS/SurveySimulation/linearJScheduler_DDPC.py +++ b/EXOSIMS/SurveySimulation/linearJScheduler_DDPC.py @@ -175,7 +175,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -244,6 +244,9 @@ def next_target(self, old_sInd, modes): # create DRM DRM = {} + # selecting appropriate koMap + koMap = self.koMaps[modes[0]['syst']['name']] + # allocate settling time + overhead time tmpCurrentTimeAbs = TK.currentTimeAbs.copy() + Obs.settlingTime + modes[0]['syst']['ohTime'] tmpCurrentTimeNorm = TK.currentTimeNorm.copy() + Obs.settlingTime + modes[0]['syst']['ohTime'] @@ -262,7 +265,7 @@ def next_target(self, old_sInd, modes): sd = None if OS.haveOcculter == True: sd = Obs.star_angularSep(TL, old_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMap,self.koTimes,modes[0]) + obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMaps,self.koTimes,modes[0]) slewTimes = Obs.calculate_slewTimes(TL, old_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -276,7 +279,7 @@ def next_target(self, old_sInd, modes): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -307,7 +310,7 @@ def next_target(self, old_sInd, modes): if len(sInds.tolist()) > 0 and Obs.checkKeepoutEnd: try: # endTimes may exist past koTimes so we have an exception to hand this case koTimeInd = np.where(np.round(endTimes[0].value)-self.koTimes.value==0)[0][0]#koTimeInd[0][0] # find indice where koTime is endTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except: sInds = np.asarray([],dtype=int) @@ -527,6 +530,9 @@ def observation_characterization(self, sInd, modes): nmodes = len(modes) + # selecting appropriate koMap + koMap = self.koMaps[modes[0]['syst']['name']] + # find indices of planets around the target pInds = np.where(SU.plan2star == sInd)[0] @@ -570,7 +576,9 @@ def observation_characterization(self, sInd, modes): startTime = TK.currentTimeAbs.copy() + mode['syst']['ohTime'] + Obs.settlingTime startTimeNorm = TK.currentTimeNorm.copy() + mode['syst']['ohTime'] + Obs.settlingTime # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + tochar[tochar] = koMap[sInd][koTimeInd] # 2/ if any planet to characterize, find the characterization times # at the detected fEZ, dMag, and WA @@ -597,8 +605,14 @@ def observation_characterization(self, sInd, modes): (endTimesNorm <= TK.OBendTimes[TK.OBnumber])) # 3/ is target still observable at the end of any char time? - if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar]) + if np.any(tochar) and Obs.checkKeepoutEnd: + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] tochars.append(tochar) intTimes_all.append(intTimes) diff --git a/EXOSIMS/SurveySimulation/linearJScheduler_DDPC_sotoSS.py b/EXOSIMS/SurveySimulation/linearJScheduler_DDPC_sotoSS.py index b89c5cc04..30d13988b 100644 --- a/EXOSIMS/SurveySimulation/linearJScheduler_DDPC_sotoSS.py +++ b/EXOSIMS/SurveySimulation/linearJScheduler_DDPC_sotoSS.py @@ -175,7 +175,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -244,6 +244,9 @@ def next_target(self, old_sInd, modes): # create DRM DRM = {} + # selecting appropriate koMap + koMap = self.koMaps[modes[0]['syst']['name']] + # allocate settling time + overhead time tmpCurrentTimeAbs = TK.currentTimeAbs.copy() + Obs.settlingTime + modes[0]['syst']['ohTime'] tmpCurrentTimeNorm = TK.currentTimeNorm.copy() + Obs.settlingTime + modes[0]['syst']['ohTime'] @@ -262,7 +265,7 @@ def next_target(self, old_sInd, modes): sd = None if OS.haveOcculter == True: sd = Obs.star_angularSep(TL, old_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMap,self.koTimes,modes[0]) + obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMaps,self.koTimes,modes[0]) slewTimes = Obs.calculate_slewTimes(TL, old_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -276,7 +279,7 @@ def next_target(self, old_sInd, modes): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -307,7 +310,7 @@ def next_target(self, old_sInd, modes): if len(sInds.tolist()) > 0 and Obs.checkKeepoutEnd: try: # endTimes may exist past koTimes so we have an exception to hand this case koTimeInd = np.where(np.round(endTimes[0].value)-self.koTimes.value==0)[0][0]#koTimeInd[0][0] # find indice where koTime is endTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except: sInds = np.asarray([],dtype=int) @@ -505,6 +508,9 @@ def observation_characterization(self, sInd, modes): nmodes = len(modes) + # selecting appropriate koMap + koMap = self.koMaps[modes[0]['syst']['name']] + # find indices of planets around the target pInds = np.where(SU.plan2star == sInd)[0] @@ -547,7 +553,9 @@ def observation_characterization(self, sInd, modes): startTime = TK.currentTimeAbs.copy() + mode['syst']['ohTime'] + Obs.settlingTime startTimeNorm = TK.currentTimeNorm.copy() + mode['syst']['ohTime'] + Obs.settlingTime # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + tochar[tochar] = koMap[sInd][koTimeInd] # 2/ if any planet to characterize, find the characterization times # at the detected fEZ, dMag, and WA @@ -570,8 +578,14 @@ def observation_characterization(self, sInd, modes): (endTimesNorm <= TK.OBendTimes[TK.OBnumber])) # 3/ is target still observable at the end of any char time? - if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar]) + if np.any(tochar) and Obs.checkKeepoutEnd: + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] tochars.append(tochar) intTimes_all.append(intTimes) diff --git a/EXOSIMS/SurveySimulation/linearJScheduler_det_only.py b/EXOSIMS/SurveySimulation/linearJScheduler_det_only.py index b3b2a36a7..676442322 100644 --- a/EXOSIMS/SurveySimulation/linearJScheduler_det_only.py +++ b/EXOSIMS/SurveySimulation/linearJScheduler_det_only.py @@ -133,7 +133,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -199,6 +199,9 @@ def next_target(self, old_sInd, mode): # create DRM DRM = {} + # selecting appropriate koMap + koMap = self.koMaps[mode['syst']['name']] + # allocate settling time + overhead time tmpCurrentTimeAbs = TK.currentTimeAbs.copy() + Obs.settlingTime + mode['syst']['ohTime'] tmpCurrentTimeNorm = TK.currentTimeNorm.copy() + Obs.settlingTime + mode['syst']['ohTime'] @@ -218,7 +221,7 @@ def next_target(self, old_sInd, mode): sd = None if OS.haveOcculter == True: sd = Obs.star_angularSep(TL, old_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMap,self.koTimes,mode) + obsTimes = Obs.calculate_observableTimes(TL,sInds,tmpCurrentTimeAbs,self.koMaps,self.koTimes,mode) slewTimes = Obs.calculate_slewTimes(TL, old_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -232,7 +235,7 @@ def next_target(self, old_sInd, mode): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -263,7 +266,7 @@ def next_target(self, old_sInd, mode): if len(sInds.tolist()) > 0 and Obs.checkKeepoutEnd: try: # endTimes may exist past koTimes so we have an exception to hand this case koTimeInd = np.where(np.round(endTimes[0].value)-self.koTimes.value==0)[0][0]#koTimeInd[0][0] # find indice where koTime is endTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except: sInds = np.asarray([],dtype=int) diff --git a/EXOSIMS/SurveySimulation/linearJScheduler_det_only_sotoSS.py b/EXOSIMS/SurveySimulation/linearJScheduler_det_only_sotoSS.py index f24e83938..f9eaa059f 100644 --- a/EXOSIMS/SurveySimulation/linearJScheduler_det_only_sotoSS.py +++ b/EXOSIMS/SurveySimulation/linearJScheduler_det_only_sotoSS.py @@ -151,7 +151,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) diff --git a/EXOSIMS/SurveySimulation/occulterJScheduler.py b/EXOSIMS/SurveySimulation/occulterJScheduler.py index d9a97ef00..570dabec5 100644 --- a/EXOSIMS/SurveySimulation/occulterJScheduler.py +++ b/EXOSIMS/SurveySimulation/occulterJScheduler.py @@ -82,7 +82,7 @@ def choose_next_target(self, old_sInd, sInds, slewTimes, intTimes): # only consider slew distance when there's an occulter if OS.haveOcculter: - angdists = Obs.star_angularSep(TL, old_sInd, sInds, dt) + angdists = [Obs.star_angularSep(TL, old_sInd, s, t) for s,t in zip(sInds,dt)] try: Obs.__getattribute__('dV_interp') @@ -115,7 +115,7 @@ def choose_next_target(self, old_sInd, sInds, slewTimes, intTimes): A_ = np.zeros((nStars,nStars)) # only consider slew distance when there's an occulter if OS.haveOcculter: - angdists_ = np.array([Obs.star_angularSep(TL, s, sInds, dt) for s in range(len(sInds))])*u.d + angdists_ = np.array([Obs.star_angularSep(TL, s, sInds, t) for s,t in zip(sInds,dt)]) dVs_= np.array([Obs.dV_interp(slewTimes[sInds[s]],angdists_[s,:]) for s in range(len(sInds))]) A_ = self.coeffs[0]*dVs_.reshape(nStars,nStars)/(0.025*Obs.dVtot.value) # add factor due to completeness diff --git a/EXOSIMS/SurveySimulation/tieredScheduler.py b/EXOSIMS/SurveySimulation/tieredScheduler.py index 77804a450..8f070bea3 100644 --- a/EXOSIMS/SurveySimulation/tieredScheduler.py +++ b/EXOSIMS/SurveySimulation/tieredScheduler.py @@ -378,7 +378,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -507,6 +507,9 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): Obs = self.Observatory TK = self.TimeKeeping + # selecting appropriate koMap + koMap = self.koMaps[char_mode['syst']['name']] + # Create DRM DRM = {} @@ -546,7 +549,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): # sd = np.zeros(TL.nStars)*u.rad # else: sd = Obs.star_angularSep(TL, old_occ_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMap, self.koTimes, char_mode) + obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMaps, self.koTimes, char_mode) slewTimes = Obs.calculate_slewTimes(TL, old_occ_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -565,7 +568,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(occ_startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds_occ_ko = occ_sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[occ_sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds_occ_ko = occ_sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[occ_sInds])[0]]# filters inds by koMap #verified against v1.35 occ_sInds = sInds_occ_ko[np.where(np.in1d(sInds_occ_ko, HIP_sInds))[0]] except:#If there are no target stars to observe sInds_occ_ko = np.asarray([],dtype=int) @@ -573,7 +576,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -640,7 +643,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): tmpIndsbool = list() for i in np.arange(len(occ_sInds)): koTimeInd = np.where(np.round(endTimes[occ_sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind occ_sInds = occ_sInds[tmpIndsbool] del tmpIndsbool except: @@ -651,7 +654,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(endTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool except: @@ -1030,6 +1033,9 @@ def observation_characterization(self, sInd, mode): SU = self.SimulatedUniverse Obs = self.Observatory TK = self.TimeKeeping + + # selecting appropriate koMap + koMap = self.koMaps[mode['syst']['name']] # find indices of planets around the target pInds = np.where(SU.plan2star == sInd)[0] @@ -1070,7 +1076,9 @@ def observation_characterization(self, sInd, mode): startTime = TK.currentTimeAbs.copy() + mode['syst']['ohTime'] + Obs.settlingTime startTimeNorm = TK.currentTimeNorm.copy() + mode['syst']['ohTime'] + Obs.settlingTime # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + tochar[tochar] = koMap[sInd][koTimeInd] # 2/ if any planet to characterize, find the characterization times if np.any(tochar): @@ -1108,7 +1116,13 @@ def observation_characterization(self, sInd, mode): # 3/ is target still observable at the end of any char time? if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar]) + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] # 4/ if yes, perform the characterization for the maximum char time if np.any(tochar): diff --git a/EXOSIMS/SurveySimulation/tieredScheduler_DD.py b/EXOSIMS/SurveySimulation/tieredScheduler_DD.py index 4b903575c..93a1e00b3 100644 --- a/EXOSIMS/SurveySimulation/tieredScheduler_DD.py +++ b/EXOSIMS/SurveySimulation/tieredScheduler_DD.py @@ -264,7 +264,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -337,6 +337,9 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # Create DRM DRM = {} + + # selecting appropriate koMap + koMap = self.koMaps[char_mode['syst']['name']] # In case of an occulter, initialize slew time factor # (add transit time and reduce starshade mass) @@ -374,7 +377,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # sd = np.zeros(TL.nStars)*u.rad # else: sd = Obs.star_angularSep(TL, old_occ_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMap, self.koTimes, char_mode) + obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMaps, self.koTimes, char_mode) slewTimes = Obs.calculate_slewTimes(TL, old_occ_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -393,7 +396,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(occ_startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds_occ_ko = occ_sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[occ_sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds_occ_ko = occ_sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[occ_sInds])[0]]# filters inds by koMap #verified against v1.35 occ_sInds = sInds_occ_ko[np.where(np.in1d(sInds_occ_ko, HIP_sInds))[0]] except:#If there are no target stars to observe sInds_occ_ko = np.asarray([],dtype=int) @@ -401,7 +404,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -465,7 +468,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): tmpIndsbool = list() for i in np.arange(len(occ_sInds)): koTimeInd = np.where(np.round(endTimes[occ_sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind occ_sInds = occ_sInds[tmpIndsbool] del tmpIndsbool except: @@ -476,7 +479,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(endTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool except: diff --git a/EXOSIMS/SurveySimulation/tieredScheduler_DD_SLSQP.py b/EXOSIMS/SurveySimulation/tieredScheduler_DD_SLSQP.py index 4958d6a55..16c1bb59a 100644 --- a/EXOSIMS/SurveySimulation/tieredScheduler_DD_SLSQP.py +++ b/EXOSIMS/SurveySimulation/tieredScheduler_DD_SLSQP.py @@ -264,7 +264,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -337,6 +337,9 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # Create DRM DRM = {} + + # selecting appropriate koMap + koMap = self.koMaps[char_mode['syst']['name']] # In case of an occulter, initialize slew time factor # (add transit time and reduce starshade mass) @@ -373,7 +376,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # sd = np.zeros(TL.nStars)*u.rad # else: sd = Obs.star_angularSep(TL, old_occ_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMap, self.koTimes, char_mode) + obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMaps, self.koTimes, char_mode) slewTimes = Obs.calculate_slewTimes(TL, old_occ_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -392,7 +395,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(occ_startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds_occ_ko = occ_sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[occ_sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds_occ_ko = occ_sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[occ_sInds])[0]]# filters inds by koMap #verified against v1.35 occ_sInds = sInds_occ_ko[np.where(np.in1d(sInds_occ_ko, HIP_sInds))[0]] except:#If there are no target stars to observe sInds_occ_ko = np.asarray([],dtype=int) @@ -400,7 +403,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -458,7 +461,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): tmpIndsbool = list() for i in np.arange(len(occ_sInds)): koTimeInd = np.where(np.round(endTimes[occ_sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind occ_sInds = occ_sInds[tmpIndsbool] del tmpIndsbool except: @@ -469,7 +472,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(endTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool except: diff --git a/EXOSIMS/SurveySimulation/tieredScheduler_DD_sotoSS.py b/EXOSIMS/SurveySimulation/tieredScheduler_DD_sotoSS.py index 290148027..498581177 100644 --- a/EXOSIMS/SurveySimulation/tieredScheduler_DD_sotoSS.py +++ b/EXOSIMS/SurveySimulation/tieredScheduler_DD_sotoSS.py @@ -260,7 +260,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -333,6 +333,9 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # Create DRM DRM = {} + + # selecting appropriate koMap + koMap = self.koMaps[char_mode['syst']['name']] # In case of an occulter, initialize slew time factor # (add transit time and reduce starshade mass) @@ -364,7 +367,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # 1 Find spacecraft orbital START positions and filter out unavailable # targets. If occulter, each target has its own START position. sd = Obs.star_angularSep(TL, old_occ_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMap, self.koTimes, char_mode) + obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMaps, self.koTimes, char_mode) slewTimes = Obs.calculate_slewTimes(TL, old_occ_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -381,7 +384,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(occ_startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds_occ_ko = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds_occ_ko = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 occ_sInds = sInds_occ_ko[np.where(np.in1d(sInds_occ_ko, HIP_sInds))[0]] except:#If there are no target stars to observe sInds_occ_ko = np.asarray([],dtype=int) @@ -389,7 +392,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -452,7 +455,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): tmpIndsbool = list() for i in np.arange(len(occ_sInds)): koTimeInd = np.where(np.round(endTimes[occ_sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind occ_sInds = occ_sInds[tmpIndsbool] del tmpIndsbool except: @@ -463,7 +466,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_modes, char_mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(endTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool except: diff --git a/EXOSIMS/SurveySimulation/tieredScheduler_SLSQP.py b/EXOSIMS/SurveySimulation/tieredScheduler_SLSQP.py index 4d33d6755..4f4e0fdd8 100644 --- a/EXOSIMS/SurveySimulation/tieredScheduler_SLSQP.py +++ b/EXOSIMS/SurveySimulation/tieredScheduler_SLSQP.py @@ -377,7 +377,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -509,6 +509,9 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): # Create DRM DRM = {} + # selecting appropriate koMap + koMap = self.koMaps[char_mode['syst']['name']] + # In case of an occulter, initialize slew time factor # (add transit time and reduce starshade mass) assert OS.haveOcculter == True @@ -544,7 +547,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): # sd = np.zeros(TL.nStars)*u.rad # else: sd = Obs.star_angularSep(TL, old_occ_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMap, self.koTimes, char_mode) + obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMaps, self.koTimes, char_mode) slewTimes = Obs.calculate_slewTimes(TL, old_occ_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -563,7 +566,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(occ_startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds_occ_ko = occ_sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[occ_sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds_occ_ko = occ_sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[occ_sInds])[0]]# filters inds by koMap #verified against v1.35 occ_sInds = sInds_occ_ko[np.where(np.in1d(sInds_occ_ko, HIP_sInds))[0]] except:#If there are no target stars to observe sInds_occ_ko = np.asarray([],dtype=int) @@ -571,7 +574,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -638,7 +641,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): tmpIndsbool = list() for i in np.arange(len(occ_sInds)): koTimeInd = np.where(np.round(endTimes[occ_sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind occ_sInds = occ_sInds[tmpIndsbool] del tmpIndsbool except: @@ -649,7 +652,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(endTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool except: @@ -1012,6 +1015,9 @@ def observation_characterization(self, sInd, mode): SU = self.SimulatedUniverse Obs = self.Observatory TK = self.TimeKeeping + + # selecting appropriate koMap + koMap = self.koMaps[mode['syst']['name']] # find indices of planets around the target pInds = np.where(SU.plan2star == sInd)[0] @@ -1052,7 +1058,9 @@ def observation_characterization(self, sInd, mode): startTime = TK.currentTimeAbs.copy() + mode['syst']['ohTime'] + Obs.settlingTime startTimeNorm = TK.currentTimeNorm.copy() + mode['syst']['ohTime'] + Obs.settlingTime # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + tochar[tochar] = koMap[sInd][koTimeInd] # 2/ if any planet to characterize, find the characterization times if np.any(tochar): @@ -1090,7 +1098,13 @@ def observation_characterization(self, sInd, mode): # 3/ is target still observable at the end of any char time? if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar]) + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] # 4/ if yes, perform the characterization for the maximum char time if np.any(tochar): diff --git a/EXOSIMS/SurveySimulation/tieredScheduler_sotoSS.py b/EXOSIMS/SurveySimulation/tieredScheduler_sotoSS.py index 0cf82ef5a..3cd771176 100644 --- a/EXOSIMS/SurveySimulation/tieredScheduler_sotoSS.py +++ b/EXOSIMS/SurveySimulation/tieredScheduler_sotoSS.py @@ -365,7 +365,7 @@ def run_sim(self): self.vprint('waitTime is not None') else: startTimes = TK.currentTimeAbs.copy() + np.zeros(TL.nStars)*u.d # Start Times of Observations - observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMap,self.koTimes,self.mode)[0] + observableTimes = Obs.calculate_observableTimes(TL,np.arange(TL.nStars),startTimes,self.koMaps,self.koTimes,self.mode)[0] #CASE 2 If There are no observable targets for the rest of the mission if((observableTimes[(TK.missionFinishAbs.copy().value*u.d > observableTimes.value*u.d)*(observableTimes.value*u.d >= TK.currentTimeAbs.copy().value*u.d)].shape[0]) == 0):#Are there any stars coming out of keepout before end of mission self.vprint('No Observable Targets for Remainder of mission at currentTimeNorm= ' + str(TK.currentTimeNorm.copy())) @@ -495,6 +495,9 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): # Create DRM DRM = {} + # selecting appropriate koMap + koMap = self.koMaps[char_mode['syst']['name']] + # In case of an occulter, initialize slew time factor # (add transit time and reduce starshade mass) assert OS.haveOcculter == True @@ -530,7 +533,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): # sd = np.zeros(TL.nStars)*u.rad # else: sd = Obs.star_angularSep(TL, old_occ_sInd, sInds, tmpCurrentTimeAbs) - obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMap, self.koTimes, char_mode) + obsTimes = Obs.calculate_observableTimes(TL, sInds, tmpCurrentTimeAbs, self.koMaps, self.koTimes, char_mode) slewTimes = Obs.calculate_slewTimes(TL, old_occ_sInd, sInds, sd, obsTimes, tmpCurrentTimeAbs) # 2.1 filter out totTimes > integration cutoff @@ -547,7 +550,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): # 2.5 Filter stars not observable at startTimes try: koTimeInd = np.where(np.round(occ_startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds_occ_ko = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds_occ_ko = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 occ_sInds = sInds_occ_ko[np.where(np.in1d(sInds_occ_ko, HIP_sInds))[0]] except:#If there are no target stars to observe sInds_occ_ko = np.asarray([],dtype=int) @@ -555,7 +558,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): try: koTimeInd = np.where(np.round(startTimes[0].value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] - sInds = sInds[np.where(np.transpose(self.koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 + sInds = sInds[np.where(np.transpose(koMap)[koTimeInd].astype(bool)[sInds])[0]]# filters inds by koMap #verified against v1.35 except:#If there are no target stars to observe sInds = np.asarray([],dtype=int) @@ -621,7 +624,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): tmpIndsbool = list() for i in np.arange(len(occ_sInds)): koTimeInd = np.where(np.round(endTimes[occ_sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[occ_sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind occ_sInds = occ_sInds[tmpIndsbool] del tmpIndsbool except: @@ -632,7 +635,7 @@ def next_target(self, old_sInd, old_occ_sInd, det_mode, char_mode): tmpIndsbool = list() for i in np.arange(len(sInds)): koTimeInd = np.where(np.round(endTimes[sInds[i]].value)-self.koTimes.value==0)[0][0] # find indice where koTime is endTime[0] - tmpIndsbool.append(self.koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind + tmpIndsbool.append(koMap[sInds[i]][koTimeInd].astype(bool)) #Is star observable at time ind sInds = sInds[tmpIndsbool] del tmpIndsbool except: @@ -991,6 +994,9 @@ def observation_characterization(self, sInd, mode): SU = self.SimulatedUniverse Obs = self.Observatory TK = self.TimeKeeping + + # selecting appropriate koMap + koMap = self.koMaps[mode['syst']['name']] # find indices of planets around the target pInds = np.where(SU.plan2star == sInd)[0] @@ -1031,7 +1037,9 @@ def observation_characterization(self, sInd, mode): startTime = TK.currentTimeAbs.copy() + mode['syst']['ohTime'] + Obs.settlingTime startTimeNorm = TK.currentTimeNorm.copy() + mode['syst']['ohTime'] + Obs.settlingTime # planets to characterize - tochar[tochar] = Obs.keepout(TL, sInd, startTime) + koTimeInd = np.where(np.round(startTime.value)-self.koTimes.value==0)[0][0] # find indice where koTime is startTime[0] + #wherever koMap is 1, the target is observable + tochar[tochar] = koMap[sInd][koTimeInd] # 2/ if any planet to characterize, find the characterization times if np.any(tochar): @@ -1065,7 +1073,13 @@ def observation_characterization(self, sInd, mode): # 3/ is target still observable at the end of any char time? if np.any(tochar) and Obs.checkKeepoutEnd: - tochar[tochar] = Obs.keepout(TL, sInd, endTimes[tochar]) + koTimeInds = np.zeros(len(endTimes.value),dtype=int) + for t,endTime in enumerate(endTimes.value): + if endTime > self.koTimes.value[-1]: + koTimeInds[t] = np.where(np.floor(endTime)-self.koTimes.value==0)[0][0] + else: + koTimeInds[t] = np.where(np.round(endTime)-self.koTimes.value==0)[0][0] # find indice where koTime is endTimes[0] + tochar[tochar] = koMap[sInd][koTimeInds] # 4/ if yes, perform the characterization for the maximum char time if np.any(tochar): diff --git a/EXOSIMS/ZodiacalLight/Stark.py b/EXOSIMS/ZodiacalLight/Stark.py index a0c8505e2..f75ba8bb5 100644 --- a/EXOSIMS/ZodiacalLight/Stark.py +++ b/EXOSIMS/ZodiacalLight/Stark.py @@ -165,11 +165,11 @@ def calcfZmax(self, sInds, Obs, TL, TK, mode, hashname): timeArray = [j*dt for j in range(1000)] #When are stars in KO regions - kogoodStart = np.zeros([len(timeArray),sInds.shape[0]])#replaced self.schedule with sInds - for i in np.arange(len(timeArray)): - kogoodStart[i,:] = Obs.keepout(TL, sInds, TK.currentTimeAbs+timeArray[i]*u.d)#replaced self.schedule with sInds - kogoodStart[i,:] = (np.zeros(kogoodStart[i,:].shape[0])+1)*kogoodStart[i,:] - kogoodStart[kogoodStart==0] = nan +# kogoodStart = np.zeros([len(timeArray),sInds.shape[0]])#replaced self.schedule with sInds +# for i in np.arange(len(timeArray)): +# kogoodStart[i,:] = Obs.keepout(TL, sInds, TK.currentTimeAbs+timeArray[i]*u.d)#replaced self.schedule with sInds +# kogoodStart[i,:] = (np.zeros(kogoodStart[i,:].shape[0])+1)*kogoodStart[i,:] +# kogoodStart[kogoodStart==0] = nan #Filter Out fZ where star is in KO region diff --git a/tests/Observatory/test_Observatory.py b/tests/Observatory/test_Observatory.py index 70f691aca..179e62a61 100644 --- a/tests/Observatory/test_Observatory.py +++ b/tests/Observatory/test_Observatory.py @@ -51,11 +51,10 @@ def test_init(self): Test of initialization and __init__. """ - req_atts = ['koAngleMin', 'koAngleMinMoon', 'koAngleMinEarth', 'koAngleMax', 'koAngleSmall', - 'ko_dtStep', 'settlingTime', 'thrust', 'slewIsp', 'scMass', 'dryMass', 'coMass', - 'occulterSep', 'skIsp', 'defburnPortion', 'checkKeepoutEnd', 'forceStaticEphem', - 'constTOF', 'occ_dtmin', 'occ_dtmax', 'maxdVpct', 'dVtot', 'dVmax', 'flowRate', - 'havejplephem'] + req_atts = ['koAngles_SolarPanel', 'ko_dtStep', 'settlingTime', 'thrust', 'slewIsp', + 'scMass', 'dryMass', 'coMass', 'occulterSep', 'skIsp', 'defburnPortion', + 'checkKeepoutEnd', 'forceStaticEphem', 'constTOF', 'occ_dtmin', 'occ_dtmax', + 'maxdVpct', 'dVtot', 'dVmax', 'flowRate', 'havejplephem'] for mod in self.allmods: with RedirectStreams(stdout=self.dev_null): @@ -114,11 +113,10 @@ def test_str(self): Test __str__ method, for full coverage and check that all modules have required attributes. """ - atts_list = ['koAngleMin', 'koAngleMinMoon', 'koAngleMinEarth', 'koAngleMax', 'koAngleSmall', - 'ko_dtStep', 'settlingTime', 'thrust', 'slewIsp', 'scMass', 'dryMass', 'coMass', - 'occulterSep', 'skIsp', 'defburnPortion', 'checkKeepoutEnd', 'forceStaticEphem', - 'constTOF', 'occ_dtmin', 'occ_dtmax', 'maxdVpct', 'dVtot', 'dVmax', 'flowRate', - 'havejplephem'] + atts_list = ['koAngles_SolarPanel', 'ko_dtStep', 'settlingTime', 'thrust', 'slewIsp', + 'scMass', 'dryMass', 'coMass', 'occulterSep', 'skIsp', 'defburnPortion', + 'checkKeepoutEnd', 'forceStaticEphem', 'constTOF', 'occ_dtmin', 'occ_dtmax', + 'maxdVpct', 'dVtot', 'dVmax', 'flowRate', 'havejplephem'] for mod in self.allmods: with RedirectStreams(stdout=self.dev_null): diff --git a/tests/OpticalSystem/test_OpticalSystem.py b/tests/OpticalSystem/test_OpticalSystem.py index 0801595c1..f17656d96 100644 --- a/tests/OpticalSystem/test_OpticalSystem.py +++ b/tests/OpticalSystem/test_OpticalSystem.py @@ -202,7 +202,8 @@ def test_str(self): """ atts_list = ['obscurFac','shapeFac','pupilDiam','intCutoff','dMag0','ref_dMag','ref_Time', - 'pupilArea','haveOcculter','IWA','OWA','WA0'] + 'pupilArea','haveOcculter','IWA','OWA','WA0','koAngles_Sun','koAngles_Earth', + 'koAngles_Moon','koAngles_Small'] for mod in self.allmods: with RedirectStreams(stdout=self.dev_null): diff --git a/tests/Prototypes/test_Observatory.py b/tests/Prototypes/test_Observatory.py index 7ebb79e1d..8f3d917d8 100644 --- a/tests/Prototypes/test_Observatory.py +++ b/tests/Prototypes/test_Observatory.py @@ -67,12 +67,16 @@ def starprop(x,y,z): return np.array([[-8.82544227, -5.93387264, 3.69977354],\ star_catalog = MockStarCatalog() t_ref = Time(2000.5, format='jyear') obs = self.fixture - # the routine under test - syst = {'occulter':False} - kogood = obs.keepout(star_catalog, np.arange(star_catalog.nStars),t_ref) + #generating koangle array with 1 subsystem + nSystems = 1 + koAngles = np.zeros([nSystems,4,2]) + koAngles[0,:,0] = 0 + koAngles[0,:,1] = 180 + kogood = obs.keepout(star_catalog, np.arange(star_catalog.nStars),t_ref,koAngles) # return value should be True self.assertTrue(np.all(kogood)) - self.assertEqual(len(kogood), star_catalog.nStars) + self.assertEqual(kogood.shape[0], nSystems) + self.assertEqual(kogood.shape[1], star_catalog.nStars) def test_cent(self): diff --git a/tests/Prototypes/test_OpticalSystem.py b/tests/Prototypes/test_OpticalSystem.py index 9c2dba9bf..df463ebd3 100644 --- a/tests/Prototypes/test_OpticalSystem.py +++ b/tests/Prototypes/test_OpticalSystem.py @@ -496,7 +496,11 @@ def outspec_compare(self, outspec1, outspec2): # which are lists of dictionaries for (d1, d2) in zip(outspec1[k], outspec2[k]): for kk in d1: - self.assertEqual(d1[kk], d2[kk]) + if kk.split("_")[0] == 'koAngles': + self.assertEqual(d1[kk][0], d2[kk][0]) + self.assertEqual(d1[kk][1], d2[kk][1]) + else: + self.assertEqual(d1[kk], d2[kk]) else: # these are strings or numbers, but not Quantity's, # because the outspec does not contain Quantity's diff --git a/tests/Prototypes/test_TargetList.py b/tests/Prototypes/test_TargetList.py index 778e9228e..46bef58fb 100644 --- a/tests/Prototypes/test_TargetList.py +++ b/tests/Prototypes/test_TargetList.py @@ -4,6 +4,7 @@ from EXOSIMS.Prototypes import TargetList import numpy as np from astropy import units as u +from astropy.time import Time from tests.TestSupport.Info import resource_path from tests.TestSupport.Utilities import RedirectStreams import json @@ -249,6 +250,53 @@ def test_completeness_specs(self): tl = TargetList.TargetList(**spec2) self.assertNotEqual(tl.PlanetPopulation.__class__.__name__,tl.Completeness.PlanetPopulation.__class__.__name__) + def test_starprop(self): + """ + Test starprop outputs + """ + + # setting up 1-dim and multi-dim arrays + timeRange = np.arange( 2000.5, 2019.5, 5 ) # 1x4 time array + timeRangeEql = np.linspace( 2000.5, 2019.5, 5 ) # 1x5 time array, same size as sInds later + + # time Quantity arrays + t_ref = Time(timeRange[0], format='jyear') # 1x1 time array + t_refArray = Time(timeRange, format='jyear') # 1x4 time array + t_refEqual = Time(timeRangeEql, format='jyear') # 1x5 time array, equal to sInds size + t_refCopy = Time(np.tile(timeRange[0],5), format='jyear') # 1x5 time array, all elements are equal + + # sInd arrays + sInd = np.array([0]) + sInds = np.array([0,1,2,3,4]) + + # testing Static Stars (set up as a default) + r_targSSBothSingle = self.targetlist.starprop(sInd ,t_ref) # should be 1x3 + r_targSSMultSinds = self.targetlist.starprop(sInds,t_ref) # should be 5x3 + r_targSSMultBoth = self.targetlist.starprop(sInds,t_refArray) # should be 5x4x3 + r_targSSEqualBoth = self.targetlist.starprop(sInds,t_refEqual) # should be 5x3 + r_targSSCopyTimes = self.targetlist.starprop(sInd ,t_refCopy) # should be 1x3 (equal defaults to 1) + + self.assertEqual( r_targSSBothSingle.shape , (1,3) ) + self.assertEqual( r_targSSMultSinds.shape , (sInds.size,3) ) + self.assertEqual( r_targSSMultBoth.shape , (t_refArray.size,sInds.size,3) ) + self.assertEqual( r_targSSEqualBoth.shape , (sInds.size,3) ) + self.assertEqual( r_targSSCopyTimes.shape , (1,3) ) + + # testing without Static Stars + self.targetlist.starprop_static = None + r_targBothSingle = self.targetlist.starprop(sInd ,t_ref) + r_targMultSinds = self.targetlist.starprop(sInds,t_ref) + r_targMultTimes = self.targetlist.starprop(sInd ,t_refArray) # should be 5x3 + r_targMultBoth = self.targetlist.starprop(sInds,t_refArray) + r_targEqualBoth = self.targetlist.starprop(sInds,t_refEqual) + r_targCopyTimes = self.targetlist.starprop(sInd ,t_refCopy) + + self.assertEqual( r_targBothSingle.shape , (1,3) ) + self.assertEqual( r_targMultSinds.shape , (sInds.size,3) ) + self.assertEqual( r_targMultTimes.shape , (t_refArray.size,3) ) + self.assertEqual( r_targMultBoth.shape , (t_refArray.size,sInds.size,3) ) + self.assertEqual( r_targEqualBoth.shape , (sInds.size,3) ) + self.assertEqual( r_targCopyTimes.shape , (1,3) ) if __name__ == '__main__':