%****************************************************************** %Create info files for Deep Impact ITS, MRI, HRI %****************************************************************** % This program creates an info file for each fit file. For example the % file N2516167681.INFO is created for the fit file N2516167681.FIT. % % This program takes the following input arguments: % % 1. kernelfiles - a SPICE meta-kernel file containing the paths to the kernel files % 2. body - IAU name of the target body, all caps % 3. scid - SPICE spacecraft id % 4. instrframe - SPICE instrument frame name % 5. fitstimekeyword - FITS header time keyword, UTC assumed % 6. input file list - path to file in which all image files are listed % 7. output folder - path to folder where infofiles should be saved to % 8. output file list - path to file in which all files for which an infofile was % created will be listed along with their start times. function create_info_files() % metakernel = ; % body = ; scid = -140; % instrframe = ; % fitstimekeyword = ; % inputfilelist = ; % infofilefolder = ; % outputfilelist = ; %SPICE kernels cspice_furnsh('C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\SPICE\kernels\sclk\DIF_SCLKSCET.00121.tsc'); cspice_furnsh('C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\SPICE\kernels\lsk\naif0012.tls'); %Sumfiles and images. Sumfile names are by MET. sumfileDir = 'C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\testData\sumfiles'; imageDir = 'C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\testData\images'; %Translate file names between MET naming and EPOXI YYMMDDHH naming %convention. itsFileTranslate = 'C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\v3_filenameByYYMMDDHH_FromCarolynEmail\ITS\dii-c-its-3_4-9p-encounter-v3.0\document\its_translate_product_id.tab.xlsx'; hriFileTranslate = 'C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\v3_filenameByYYMMDDHH_FromCarolynEmail\HRI\dif-c-hriv-3_4-9p-encounter-v3.0\document\hriv_translate_product_id.tab'; mriFileTranslate = ''; %Image file index. This tells us the path to the EPOXI-named file. itsIndex = 'C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\v3_filenameByYYMMDDHH_FromCarolynEmail\ITS\dii-c-its-3_4-9p-encounter-v3.0\index\index.tab.xlsx'; hriIndex = 'C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\v3_filenameByYYMMDDHH_FromCarolynEmail\HRI\dif-c-hriv-3_4-9p-encounter-v3.0\index\index.tab'; mriIndex = 'C:\Users\nguyel1\Projects\SBMT\CarolynErnstGrant\data\DeepImpact\v3_filenameByYYMMDDHH_FromCarolynEmail\MRI\dif-c-hriv-3_4-9p-encounter-v3.0\index\index.tab'; mriImages = getImages(scid, sumfileDir, 'IV', itsFileTranslate, itsIndex) cspice_furnsh(metakernel); fitfiles = textscan(inputfilelist, '%s'); % Image list % ofstream fout(outputfilelist.c_str()); % if (!fout.is_open()) % cerr << "Error: Unable to open file for writing" << endl; % exit(1); % end % % for i=0:numel(fitfiles) % % getEt(fitfiles[i], sclkkey, utc, et, scid); % % getScPositionAndOrientation(et, body, scid, instr, scposb, boredir, updir, frustum); % % getSunPosition(et, body, sunPosition); % % const size_t last_slash_idx = fitfiles[i].find_last_of("\\/"); % if (std::string::npos != last_slash_idx) % { % fitfiles[i].erase(0, last_slash_idx + 1); % } % int length = fitfiles[i].size(); % string infofilename = infofilefolder + "/" % + fitfiles[i].substr(0, length-4) + '.INFO'; % saveInfoFile(infofilename, utc, scposb, boredir, updir, frustum, sunPosition); % cout << "created " << infofilename << endl; % % fout << fitfiles[i].substr(0, length) << " " << utc << endl; % % end % cout << "done." << endl; end %Find the FIT files associated with the input sumfiles. % %Read in all the sumfiles in the given directory and parse out the %image MET from the filename. Convert that to UTC YYMMDDHH and find %the matching image in imageDir. Add that image to fitfiles. function fitfiles = getImages(scid, sumfilesDir, imagerFilenamePrefix, filenameTranslator, indexFile) sumfiles = dir(sumfilesDir); for k = 1:length(sumfiles) filename = sumfiles(k).name; [pathstr,name,ext] = fileparts(filename); %Is the current file a sumfile? if strcmp(ext, '.SUM') || strcmp(ext, '.sum') sprintf('current sumfile is %s', name) found = strfind(name, imagerFilenamePrefix); %Is it for the requested DeepImpact imager? if ~isempty(found) %Parse out the MET from the sumfile name token = strsplit(name, imagerFilenamePrefix); MET = token(2); %Find the image in the translator file % ref = dlmread(filenameTranslator); t = readtable(filenameTranslator,'ReadVariableNames',false); row = find(strcmp(t, name, 1, 'first'); metNames = t(:,3); fid = fopen(filenameTranslator, 'rt'); FormatString = repmat('%s',1,5); % read the entire file, if not too big s = textscan(fid, '%s', 'delimiter', '\n'); % search for your Region: idx1 = find(strcmp(s{1}, MET), 1, 'first'); idx2 = find(strcmp(s{1}, MET), 1, 'last'); s{1}(idx1+1:idx2-1) % and read from s{1}(idx1+1:idx2-1) the values using textscan again ... disp(s{1}(idx1)) fclose(fid); % %Convert MET to UTC YYMMDDHH % et = cspice_scs2e( scid, MET ); % % et2utc_c(et, 'ISOC', 3, 25, utc); % output = cspice_timout( et, 'YYYY-MM:DD:HR:MN:SC ::UTC' ) % output = cspice_timout( et, 'YYYYMMDDHR ::UTC' ) end end end end % void splitFitsHeaderLineIntoKeyAndValue(const string& line, string& key, % string& value) { % key = line.substr(0, 8); % trim(key); % value = line.substr(10); % size_t found = value.find_last_of("/"); % if (found != string::npos) % value = value.substr(0, found); % trim(value); % removeSurroundingQuotes(value); % trim(value); % } % % void getFieldsFromFitsHeader(const string& labelfilename, string& startmet, % string& stopmet, string& target, int& naxis1, int& naxis2) % { % ifstream fin(labelfilename.c_str()); % % if (fin.is_open()) { % char buffer[81]; % string str; % string key; % string value; % % for (int i = 0; i < 100; ++i) { % fin.read(buffer, 80); % buffer[80] = '\0'; % str = buffer; % splitFitsHeaderLineIntoKeyAndValue(str, key, value); % % if (key == "NAXIS1") { % naxis1 = atoi(value.c_str()); % } else if (key == "NAXIS2") { % naxis2 = atoi(value.c_str()); % } else if (key == "SPCSCLK") { % startmet = value; % stopmet = value; % } else if (key == "TARGET") { % target = value; % } % } % } else { % cerr << "Error: Unable to open file '" << labelfilename << "'" << endl; % exit(1); % } % % fin.close(); % } % % void getStringFieldFromFitsHeader(const string& fitfile, % string key, % string& value ) % { % ifstream fin(fitfile.c_str()); % % if (fin.is_open()) % { % string str; % char buffer[81]; % string currkey; % string val; % while(true) % { % fin.read(buffer, 80); % buffer[80] = '\0'; % str = buffer; % splitFitsHeaderLineIntoKeyAndValue(str, currkey, val); % if (currkey == key) % { % value = val.c_str(); % break; % } % } % } % else % { % cerr << "Error: Unable to open file '" << fitfile << "'" << endl; % exit(1); % } % % fin.close(); % } % % void getEt(const string& fitfile, % string sclkkey, % string& utc, % double& et, % const char* scid) % { % int instid; % int found; % bodn2c_c(scid, &instid, &found); % if (failed_c() || !found) % return; % % ifstream fin(fitfile.c_str()); % % if (fin.is_open()) % { % string str; % int len = 25; % char utcchar[len]; % getStringFieldFromFitsHeader(fitfile, sclkkey, str); % % scs2e_c(instid, str.c_str(), &et); % et2utc_c(et, "ISOC", 3, len, utcchar); % utc = utcchar; % } % else % { % cerr << "Error: Unable to open file '" << fitfile << "'" << endl; % exit(1); % } % % fin.close(); % } % % /* % This function calculates spacecraft position and orientation. % % Input: % et: Ephemeris time % body: Name of celestial body which is the target % obs: Name of observing spacecraft % instrFrame: NAIF frame ID of instrument on the observing spacecraft % % Output: % scposbf: Spacecraft position in bodyframe coordinates % boredir: Boresight direction in bodyframe coordinates % updir: % frustum: Field of view boundary corner vectors in bodyframe coordinates % % */ % void getScPositionAndOrientation(double et, string body, const char* obs, const char* instrFrame, % double scposbf[3], double boredir[3], double updir[3], double frustum[12]) % { % double lt; % double targpos[3]; %sc to target vector in j2000 % double scpos[3]; %target to sc vector in j2000 % double inst2inert[3][3], inert2bf[3][3], inst2bf[3][3], rot[3][3]; % char shape[32]; % char frame[32]; % double bsight [3]; % int n; % double bounds [MAXBND][3]; % double boundssbmt [MAXBND][3]; % % The celestial body is the target when dealing with light time % const char* target = body.c_str(); % string ref = string("IAU_") + body.c_str(); % const char* abcorr = "LT+S"; % const char* inertframe = "J2000"; % SpiceInt instid; % double tmpvec[3]; % % /* % * Compute the apparent position of the center of the target body % * as seen from the spacecraft at the epoch of observation (et), % * and the one-way light time from the target to the spacecraft. % */ % spkpos_c(target, et, inertframe, abcorr, obs, targpos, <); % if (failed_c()) { % cerr << "Failed spkpos" << endl; % return; % } % % /* % * Get the position of the observer. This is just the negative of the % * spacecraft-target vector using vminus(). Note that this is _NOT_ % * the same as the apparent position of the spacecraft as seen from % * the target! % */ % vminus_c(targpos, scpos); % % /* % * Get the coordinate transformation from instrument to % * inertial frame at time ET % */ % pxform_c(instrFrame, inertframe, et, inst2inert); % if (failed_c()) { % cout << "Failed pxform1" << endl; % return; % } % % /* % * Get field of view boresight and boundary corners % */ % namfrm_c(instrFrame, &instid); % if (failed_c()) { % cout << "Failed namfrm" << endl; % return; % } % getfov_c(instid, MAXBND, WDSIZE, WDSIZE, shape, frame, bsight, &n, bounds); % if (failed_c()) { % cout << "Failed getfov" << endl; % return; % } % % /* % * There is a 180 degree rotation about the boresight in the % * data, not present in the sumfiles (so can't correct in sbmt code). % */ % axisar_c(bsight, pi_c(), rot); % mxm_c(inst2inert, rot, inst2inert); % % /* % * Get the coordinate transformation from inertial to % * body-fixed coordinates at ET minus one light time ago. % */ % pxform_c(inertframe, ref.c_str(), et - lt, inert2bf); % if (failed_c()) { % cout << "Failed pxform2" << endl; % return; % } % % /* % * transform scpos vector into body-fixed from j2000 frame % */ % mxv_c(inert2bf, scpos, scposbf); % % /* % * Compute complete transformation to go from % * instrument-fixed coordinates to body-fixed coords % */ % mxm_c(inert2bf, inst2inert, inst2bf); % % %swap the boundary corner vectors so they are in the correct order for SBMT % %getfov returns them counter-clockwise starting in the +X,+Y quadrant. % %SBMT expects them in the following order (quadrants): -X,-Z -> +X,-Z -> -X,+Z -> +X,+Z % %So the vector index mapping is % %SBMT SPICE % % 0 1 % % 1 0 % % 2 2 % % 3 3 % boundssbmt[0][0] = bounds[1][0]; % boundssbmt[0][1] = bounds[1][1]; % boundssbmt[0][2] = bounds[1][2]; % boundssbmt[1][0] = bounds[0][0]; % boundssbmt[1][1] = bounds[0][1]; % boundssbmt[1][2] = bounds[0][2]; % boundssbmt[2][0] = bounds[2][0]; % boundssbmt[2][1] = bounds[2][1]; % boundssbmt[2][2] = bounds[2][2]; % boundssbmt[3][0] = bounds[3][0]; % boundssbmt[3][1] = bounds[3][1]; % boundssbmt[3][2] = bounds[3][2]; % % %transform boresight into body frame. % mxv_c(inst2bf, bsight, boredir); % % %transform boundary corners into body frame and pack into frustum array. % int k = 0; % for (int i=0; i % #include % #include % #include "SpiceUsr.h" % % using namespace std; % % #define TIMFMT "YYYY-MON-DD HR:MN:SC.###::UTC (UTC)" % #define TIMLEN 41 % % /* % This function computes the position of the target in the observer body frame, at the time the spacecraft observed the body. % % Input: % et: Ephemeris time when an image of the body was taken % spacecraft: Name of the spacecraft that took the image % observerBody: Name of observer body (e.g. EROS, PLUTO, PHOEBE) % targetBody: Name of target body (e.g. SUN, EARTH). SPICE names can be found at % https://naif.jpl.nasa.gov/pub/naif/toolkit_docs/FORTRAN/req/naif_ids.html % % Output: % bodyToSc: The position of the target in observer body-fixed coordinates corrected for light time % velocity: The velocity of the target in observer body-fixed coordinates corrected for light time % */ % void getTargetState(double et, const char* spacecraft, const char* observerBody, const char* targetBody, double bodyToTarget[3], double velocity[3]) % { % double lt, notUsed[6], bodyToTargetState[6]; % const char* abcorr = "LT+S"; % string bodyFrame = string("IAU_") + observerBody; % % /* % * Compute the apparent state of the center of the observer body % * as seen from the spacecraft at the epoch of observation (et), % * and the one-way light time from the observer to the spacecraft. % * Only the returned light time will be used from this call, as % * such, the reference frame does not matter here. Use the body % * fixed frame. % */ % spkezr_c(observerBody, et, bodyFrame.c_str(), abcorr, spacecraft, notUsed, <); % if (failed_c()) { % cerr << "Failed spkezr" << endl; % return; % } % % /* % * Back up the time at the observer body by the light time to the % * spacecraft. This is the time that light from the target body was % * received at the observer body when the spacecraft took the image. % * It is the time at the observer body. Now simply get the position % * of the target at this time, as seen from the observer body, in % * the observer body frame. % */ % spkezr_c(targetBody, et - lt, bodyFrame.c_str(), abcorr, observerBody, bodyToTargetState, <); % if (failed_c()) { % cerr << "Failed spkezr" << endl; % return; % } % % /* % * Assign the output variables. % */ % bodyToTarget[0] = bodyToTargetState[0]; % bodyToTarget[1] = bodyToTargetState[1]; % bodyToTarget[2] = bodyToTargetState[2]; % velocity[0] = bodyToTargetState[3]; % velocity[1] = bodyToTargetState[4]; % velocity[2] = bodyToTargetState[5]; % } % % % % % #include % #include % #include % #include % #include % #include "SpiceUsr.h" % % using namespace std; % % #define TIMFMT "YYYY-MON-DD HR:MN:SC.###::UTC (UTC)" % #define TIMLEN 41 % % /* % This function computes the position of the spacecraft in the observer body frame, at the time the spacecraft observed the body. % % Input: % et: Ephemeris time when the image was taken % observerBody: Name of observer body (e.g. EROS, PLUTO) % spacecraft: NAIF SPICE name of spacecraft (e.g. NEAR, NH). These can be found at % https://naif.jpl.nasa.gov/pub/naif/toolkit_docs/FORTRAN/req/naif_ids.html % % Output: % bodyToSc: The position of the spacecraft in observer body-fixed coordinates corrected for light time % velocity: The velocity of the spacecraft in observer body-fixed coordinates corrected for light time % */ % void getSpacecraftState(double et, const char* observerBody, const char* spacecraft, double bodyToSc[3], double velocity[3]) % { % double lt, scToBodyState[6]; % const char* abcorr = "LT+S"; % string bodyFrame = string("IAU_") + observerBody; % % /* % * Compute the apparent state of the body as seen from the % * spacecraft at the epoch of observation, in the body-fixed % * frame, corrected for stellar aberration and light time. % * Note that the time entered is the time at the spacecraft, % * who is the observer. % */ % spkezr_c(observerBody, et, bodyFrame.c_str(), abcorr, spacecraft, scToBodyState, <); % if (failed_c()) { % cerr << "Failed spkezr" << endl; % return; % } % % /* % * The state of the spacecraft (apparent position and velocity) % * relative to the body is just the negative of this state. Note % * that this is not the same as the apparent position and velocity % * of the spacecraft as seen from the body at time et, because et % * is the time at the spacecraft not the body. % */ % bodyToSc[0] = -scToBodyState[0]; % bodyToSc[1] = -scToBodyState[1]; % bodyToSc[2] = -scToBodyState[2]; % velocity[0] = -scToBodyState[3]; % velocity[1] = -scToBodyState[4]; % velocity[2] = -scToBodyState[5]; % } % % % % % /*--------------------------------------*/ % /* The light time calcuation derived by */ % /* Howard Taylor for SBMT NH LORRI data */ % /* works only for the spacecraft */ % /* position. Corrected functions for */ % /* non-spacecraft bodies are below. */ % /*--------------------------------------*/ % % getSpacecraftState(et, body, sc, scposb, scvelb); % getTargetState(et, sc, body, "SUN", sunpos, unused); % getTargetState(et, sc, body, "EARTH", earthpos, unused); % % /*--------------------------------------*/ % /* Both of these calculations correctly */ % /* convert et to utc. The former outputs*/ % /* time in YYY-MM-DD format, the latter */ % /* in day of year format. */ % /*--------------------------------------*/ % % et2utc_c(et, "ISOC", 3, 25, utc); % timout_c ( et, "YYYY-DOYTHR:MN:SC.### ::UTC", 32, utc );