CAMP FCW trial analysis, updated for VIRTTEX implementation
The VIRTTEX analysis is based on the program written for the NADS data by me in 2004. However, it is different in three key respects:
1) The VIRTTEX sessions contain up to 34 trial runs per data file. The NADS data had one trial run per data file. This program uses the facilities in the stateflag classes to segment the data.
2) The VIRTTEX trials contain some additional factors related to visual resolution. These are recorded in the output files for future ANOVA's
3) The code and style have been updated to the current version of the timeSeries class library
-- Jeff Greenberg May, 2005
Copyright Jeff Greenberg.
This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation; either version 2.1 of the License, or (at your option) any later version.
This library is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public License along with this library; if not, write to the Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation; either version 2.1 of the License, or (at your option) any later version.
This library is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public License along with this library; if not, write to the Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
#include <iostream> #include <fstream> #include <list> #include <string> #include <vector> #include <limits> using namespace std; #include "dataFrame.h" #include "dataFile.h" #include "timeSeries.h" #include "timeSeriesArray.h" #include "lookupTable.h" #include "TTC2.h" #include "filter.h" #include "measures.h" #include "stateFlag.h" #include "camp_utilities.h" using namespace std; const double MPHPERFPS=0.68182; const double FTPERMETER=3.28084; const double NoCrash = 999.99; std::vector<measures> process_file(dataFile& data, lookupTable& t, string& s); iir_t *init_lowpass_20(); iir_t *init_lowpass_50(); iir_t *init_cheby1_20(); iir_t *init_cheby1_50(); dataFile openDataFile(std::string fileName,int filenum,int argc); template <class T> void filtfilt(iir_t *lpf,timeSeries<T>& POVspeed_filt); template <class T> int findBrakingReactionFrame(timeSeries<T> SVaccel, timeSeries<T> SVspeed, timeSeries<T> throttle, timeSeries<T> brake); template <class T> int findSteeringReactionFrame(timeSeries<T> SVaccel, timeSeries<T> swavel, timeSeries<T> swa, timeSeries<T> lanePos, timeSeries<T> SVYawRate); void outputMeasures(const std::string& measuresFile, const std::string& headerFile, const std::vector<measures>& dependentMeasures, const int filenum); vector<state> findAllLevels(stateFlag& sf); vector<state> findAllScenarios(stateFlag& sf); vector<state> validScenarios(vector<state>& s, timeSeries<float>& subno, stateFlag& kin, timeSeries<int>& mode, timeSeries<int>& effort, timeSeries<int>& res, timeSeries<int>& aa); bool sameScenario(state& s1, state& s2, stateFlag& kinematic_cond, timeSeries<int>& mode, timeSeries<int>& effort, timeSeries<int>& resolution, timeSeries<int>& aa); string computeSubjectNo(timeSeries<float> subj); string intToCondition(const int c); string intToResolution(const int c); string intToAA(const int c); float sampleTime(string s); #ifdef WIN32 string convertLink(const char *fname); int analysisMain (int argc, char * argv[]) #else int main (int argc, char * argv[]) #endif { // get files from the command line vector<string> fileList = parseFileList(argc, argv); int nFiles = fileList.size(); // define a lookup table for subject meta-data string subjectTablePath="subject_data.txt"; try { static lookupTable subjectTable(subjectTablePath,"Subject"); for (int filenum=0;filenum<nFiles;++filenum) { try { string dataFileName = fileList[filenum]; dataFile data = openDataFile(dataFileName, filenum+1, nFiles); vector<measures> dependentMeasures=process_file(data,subjectTable, dataFileName); outputMeasures(dataFileName, "camp_virttex.tit", // label file name dependentMeasures, filenum+1); } // these are all errors that cause us to skip a data file // but not abort the run catch (tsException& tse) { tse.debug_print(); cerr << "skipping to next file..." << endl; } catch (dataFile::dataInputExceptions& die) { die.debug_print(); cerr << "skipping to next file..." << endl; } catch (lookupTable::lookupExceptions& lookupError) { lookupError.debug_print(); cerr << "skipping to next file..." << endl; } } } catch (lookupTable::lookupExceptions& lookupError) { // if we can't read a lookup file we have to just exit lookupError.debug_print(); exit(1); } } dataFile openDataFile(string fileName,int nFile,int maxFiles) { //const double Ts=0.05; #ifdef WIN32 fileName=convertLink(fileName.c_str()); #endif dataFrame expDataFrame(removePathExt(fileName)+".hdr"); string dataFileName(removePathExt(fileName)+".bin"); cout << endl; cout << "reading file " << dataFileName << " " << nFile << " of "; cout << maxFiles << "..." << flush; //open the header file and get the sample time float Ts = sampleTime(removePathExt(fileName)+".hdr"); // load the entire data file into memory dataFile data(fileName,expDataFrame,Ts,dataFile::BigEndian); cout << "done" << endl; return data; } vector<measures> process_file(dataFile& data, lookupTable& subjectTable, string& filename) { // get the variables we're interested in timeSeries<float> tapa("TAPA",data); // remove the closed plate angle from the throttle position float minTap = tapa[tapa.min()]; timeSeries<float> throttle = tapa - minTap; // create an aggregate braking torque timeSeries<float> b1("BRKTQ(1)",data); timeSeries<float> b2("BRKTQ(2)",data); timeSeries<float> b3("BRKTQ(3)",data); timeSeries<float> b4("BRKTQ(4)",data); timeSeries<float> brake = (b1+b2+b3+b4)/4.0; timeSeries<float> steerin("STEERIN",data); // VIRTTEX has a -2.5deg bias in the steering // correct it out here //timeSeries<float> swa = steerin - steerin.mean(); timeSeries<float> swa = steerin; timeSeries<float> pov_xpos("POV_XPOS",data); // pov position in feet timeSeries<float> sv_xpos("VEHYPOS",data); // sv posiition in meters const float bumper_offset=15.2; //avg of 3 test runs in VIRTTEX timeSeries<float> range = pov_xpos - sv_xpos*FTPERMETER - bumper_offset; timeSeries<float> SVYawRate("VEHOMEGAZ",data); timeSeries<float> lanepos("LANEPOS",data); // calculate steering wheel velocity timeSeries<float> swavel = swa.differentiate(); timeSeries<float> swavel_f = swavel; int sampleRate = int(floor((1.0/swavel.sampleTime()) + 0.5)); iir_t *lpf; if (sampleRate == 20) { lpf = init_cheby1_20(); } else { lpf = init_cheby1_50(); } filtfilt(lpf, swavel_f); free_iir_filter(lpf); // free the filter // the SVaccel is in ft/s, since this is (mostly) what NADS uses. We'll keep to ft // and s for the kinematics here to be consistent timeSeries<float> accelx("VEHRPACCELX",data); timeSeries<float> SVaccel = accelx*FTPERMETER; timeSeries<float> vmph("VMPH",data); timeSeries<float> SVspeed = vmph/MPHPERFPS; timeSeries<float> POVspeed("POV_SPEED",data); // POV accel is not directly available, so once again, we calculate it. // Use a 10Hz, acausal, phaseless filter to smooth POVspeed first, since the NADS // SCC objects are updated at a much lower sample frequency timeSeries<float>POVspeed_filt = POVspeed; if (sampleRate == 20) { lpf = init_lowpass_20(); } else { lpf = init_lowpass_50(); } filtfilt(lpf,POVspeed_filt); timeSeries<float>POVaccel = POVspeed_filt.differentiate(); free_iir_filter(lpf); // extract the lateral lane position timeSeries<float>SVLatLanePos = lanepos*FTPERMETER; int reactionFrame; timeSeries<float> accely("VEHRPACCELY",data); timeSeries<float> SVLatAccel = accely*FTPERMETER; stateFlag scenarioFlag = getVirttexStateFlag("SCENARIO_NUMBER",data); timeSeries<float> effort_float("EFFORT",data); //timeSeries<int> effort = timeSeries<int>(effort_float); timeSeries<int> effort = effort_float; timeSeries<float> scenario_type_float("SCENARIO_TYPE",data); timeSeries<int> scenario_type = timeSeries<int>(scenario_type_float); timeSeries<float> resolution_float("RESOLUTION",data); timeSeries<int> resolution = timeSeries<int>(resolution_float); timeSeries<float> aa_float("AA",data); timeSeries<int> aa = timeSeries<int>(aa_float); timeSeries<float> subno("SUBNO",data); string subject = computeSubjectNo(subno); cout << "subject " << subject << endl; // iterate over the scenarios in the experiment and compute the relavent // measures string subjectType = subjectTable.find<string>(subject,"Mode"); cout << "test mode: " << subjectType << endl; string gender = subjectTable.find<string>(subject,"Gender"); string ageGroup = subjectTable.find<string>(subject,"Age"); cout << "gender = " << gender << endl; cout << "age = " << ageGroup << endl; string motionType = subjectTable.find<string>(subject,"Motion"); vector<measures> event_measures; vector<state> scenarios; // if a '.scn' file exists read the scenario state definitions from there // otherwise figure them out from the data file string scnFileName = removePathExt(filename)+".scn"; ifstream scnFile(scnFileName.c_str(),ios::in); if (scnFile.is_open()) { state s(0,0); while (!scnFile.eof()) { scnFile >> s; scenarios.push_back(s); }; } else { vector<state> allScenarios = findAllScenarios(scenarioFlag); // prune out duplicates and test runs before beginning the analyis scenarios = validScenarios(allScenarios,subno,scenarioFlag,scenario_type, effort,resolution,aa); } if (getPrintState()) { cout << "============================" << endl; for (int i=0;i<scenarios.size();i++) { cout << scenarios[i] << endl; } cout << "============================" << endl; } for (int ns=0;ns<scenarios.size();ns++) { int endIndex = scenarios[ns].end(); bool brakingTrial = (scenario_type[endIndex] == 0); if (brakingTrial) { reactionFrame = findBrakingReactionFrame(SVaccel.range(scenarios[ns]), SVspeed.range(scenarios[ns]), throttle.range(scenarios[ns]), brake.range(scenarios[ns])); } else { reactionFrame = findSteeringReactionFrame(SVLatAccel.range(scenarios[ns]), swavel.range(scenarios[ns]), swa.range(scenarios[ns]), SVLatLanePos.range(scenarios[ns]), SVYawRate.range(scenarios[ns])); } cout << (brakingTrial ? "braking " : "steering ") << "onset frame= " ; cout << reactionFrame << endl; int scenario_number = scenarioFlag[endIndex]; bool stationaryTrial = (scenario_number >= 15) && (scenario_number <= 17); bool decelTrial = (scenario_number >= 1) && (scenario_number <= 9); bool constantDeltaV = (scenario_number >= 10) && (scenario_number <= 14); string instruction = (effort[endIndex] ? "hard" : "normal"); // Calculates TTC2 double ttc; try { ttc = TTC2(SVspeed[reactionFrame], POVspeed[reactionFrame], 0.0, //SVaccel (decelTrial ? POVaccel[reactionFrame] : 0.0), range[reactionFrame]); } catch (TTC::TTCExceptions& ttce) { ttc = measures::MISSING; } if (ttc == NoCrash) ttc = measures::MISSING; // Looking for mininum TTC2 from onset reactionFrame to end of run int i, finalFrame, finalFrameLC; double tmp2, minttc2; double min_seen2 = numeric_limits<double>::max(); if(brakingTrial) { finalFrame = scenarios[ns].end(); } else { timeSeries<float> lp = lanepos.range(reactionFrame,scenarios[ns].end()); const float minLC = 1.5; // lane pos to just clear the POV finalFrameLC = find(lp<= minLC).front(); finalFrame = finalFrameLC; } for (i = reactionFrame; i < finalFrame; ++i){ try { tmp2 = TTC2(SVspeed[i], //calls TTC2 function in TTC2.cpp program POVspeed[i], 0.0, //SVaccel (decelTrial ? POVaccel[reactionFrame] : 0.0), range[i]); } catch (TTC::TTCExceptions& ttce) { break; //if an exception error is found break out of for loop } if (tmp2 < min_seen2) { min_seen2 = tmp2; } } // iterating for minimum TTC2 here minttc2 = min_seen2; if (minttc2 == NoCrash) minttc2 = measures::MISSING; // calculates TTC1 and InvTTC1 double ttc1, inv_ttc1; ttc1 = range[reactionFrame]/(SVspeed[reactionFrame] - POVspeed[reactionFrame]); inv_ttc1 = 1/ttc1; // Looking for mininum TTC1 from onset reactionFrame to end of run double tmp1, minttc1; double min_seen1 = numeric_limits<double>::max(); if(brakingTrial) { finalFrame = SVspeed.getMaxIndex(); } else { finalFrame = finalFrameLC; } for (i = reactionFrame; i < finalFrame; ++i){ try { tmp1 = TTC2(SVspeed[i], //calls TTC2 function in TTC2.cpp program POVspeed[i], 0.0, //SVAccel 0.0, //POVAccel range[i]); } catch (TTC::TTCExceptions& ttce) { break; //if an exception error is found break out of for loop } if (tmp1 < min_seen1) { min_seen1 = tmp1; } } // iterating for minimum TTC1 here minttc1 = min_seen1; if (minttc1 == NoCrash) minttc1 = measures::MISSING; // Calculates Required Deceleration const double SystemDelay=0.0; double reqDec; try { reqDec= requiredDecel(SystemDelay, SVspeed[reactionFrame], POVspeed[reactionFrame], 0.0, //SVaccel (decelTrial ? POVaccel[reactionFrame] : 0.0), range[reactionFrame])/32.174; // converted to "g's" } catch (TTC::TTCExceptions& ttce) { reqDec = measures::MISSING; } // Initial Headway @ POV Deceleration = -0.1 g int priorToOnsetFrame; vector<int> onsetIndx; double Treact, initHdwyRange=0.0, initHdwyTime = 0.0; double initSVHdwySpeed, initPOVHdwySpeed; if (decelTrial) { Treact = POVaccel.time(reactionFrame); timeSeries<float> POVbeforeOnset = POVaccel.range(Treact-7.0,Treact); onsetIndx = find(POVbeforeOnset < -0.1*32.174); if (onsetIndx.size() == 0) { cout << "POV is flat over range & SV onset occurs before POV onset."; cout << endl; continue; } else { priorToOnsetFrame = onsetIndx.front(); initHdwyRange = range[priorToOnsetFrame]; initHdwyTime = range[priorToOnsetFrame]/SVspeed[priorToOnsetFrame]; } } // Calculates Peak Deceleration during braking/steering maneuver // from reactFrame to end of maneuver double pkDecel=0.0; int pkDecelFrame; if(brakingTrial) { // finds maximum longitudinal braking acceleration finalFrame = scenarios[ns].end(); pkDecelFrame = SVaccel.range(reactionFrame,finalFrame).min(); pkDecel = SVaccel[pkDecelFrame]/32.174; } else { // finds maximum lateral steering deceleration finalFrame = finalFrameLC; pkDecelFrame = SVLatAccel.range(reactionFrame,finalFrame).min(); pkDecel = SVLatAccel[pkDecelFrame]/32.174; } // find minimum range double minRange = numeric_limits<double>::max(); if (brakingTrial) { finalFrame = scenarios[ns].end(); } else { finalFrame = finalFrameLC; } for (i=reactionFrame;i<finalFrame;i++) { if (range[i] < minRange) { minRange = range[i]; } } string condition = intToCondition(scenario_number); string res = intToResolution(resolution[scenarios[ns].end()]); string antialias = intToAA(aa[scenarios[ns].end()]); cout << condition << endl; cout << instruction << " " << res << " " << antialias << endl; cout << "TTC1 = " << ttc1 << " InvTTC1 = " << inv_ttc1 << endl; cout << "minTTC1= " << minttc1 << " minTTC2= " << minttc2 << endl; cout << "TTC2= " << ttc << endl; cout << "Initial headway = " << initHdwyRange << " ft"; cout << " Peak Dec= " << pkDecel << " g's" << endl; cout << "rDec=" << reqDec << " g's" << endl; cout << endl << flush; // create the measures measures meas; meas.entry("subno", subject); meas.entry("subtype",subjectType); meas.entry("filename",getLeafName(filename)); meas.entry("gender",gender); meas.entry("age",ageGroup); meas.entry("condition", condition); meas.entry("instruction", instruction); meas.entry("mode",(brakingTrial ? "brake" : "steer")); string kinType; if (stationaryTrial) { kinType = "stationary"; } else if (decelTrial) { kinType = "decel"; } else { kinType = "constdv"; } meas.entry("kintype",kinType); meas.entry("resolution",res); meas.entry("antialias",antialias); meas.entry("motion",motionType); meas.entry("reactionTime",SVspeed.time(reactionFrame)); meas.entry("reqdec",reqDec); meas.entry("TTC2",ttc); meas.entry("SVspeed",SVspeed[reactionFrame]); meas.entry("POVspeed",POVspeed[reactionFrame]); meas.entry("SVaccel",SVaccel[reactionFrame]); meas.entry("POVaccel",(decelTrial ? POVaccel[reactionFrame] : 0.0)); meas.entry("Range", range[reactionFrame]); meas.entry("minRange",minRange); meas.entry("minTTC2", minttc2); float SVspdinit; float POVspdinit; if (decelTrial) { if (onsetIndx.size() == 0) { initHdwyRange = measures::MISSING; initHdwyTime = measures::MISSING; SVspdinit = SVspeed[reactionFrame]*MPHPERFPS; POVspdinit = POVspeed[reactionFrame]*MPHPERFPS; } else { SVspdinit = SVspeed[priorToOnsetFrame]*MPHPERFPS; POVspdinit = POVspeed[priorToOnsetFrame]*MPHPERFPS; } } else { initHdwyRange = measures::MISSING; initHdwyTime = measures::MISSING; SVspdinit = measures::MISSING; POVspdinit = measures::MISSING; } meas.entry("initHdwyRange", initHdwyRange); meas.entry("initHdwyTime", initHdwyTime); meas.entry("SVspdinit", SVspdinit); meas.entry("POVspdinit", POVspdinit); meas.entry("peakDec",pkDecel); if (ttc1 < 0) { // negative TTC1 cases are undefined ttc1 = measures::MISSING; inv_ttc1 = measures::MISSING; minttc1 = measures::MISSING; } meas.entry("TTC1",ttc1); meas.entry("InvTTC1",inv_ttc1); meas.entry("minTTC1",minttc1); // save this measures object event_measures.push_back(meas); } return event_measures; } /* lowpass_20 is a 3rd order butterworth filter with a break frequency at 9 Hz, sampling freq. at 20 Hz */ iir_t *init_lowpass_20() { double coeff[][6]= { {1, 1.00000959440837, 0, 1, 0.726542528005365, 0}, {1, 1.99999040559162, 0.999990405683675, 1, 1.64755221570399, 0.732338917272799} }; double gain=0.729440722639082; iir_t *f; f=init_iir_filter(coeff, gain, 2); return f; } /* lowpass_50 is a 3rd order butterworth filter with a break frequency at 10 Hz, sampling freq. at 50 Hz */ iir_t *init_lowpass_50() { double coeff[][6]= { {1, 1.00000959440837, 0, 1, -0.158384440324536, 0}, {1, 1.99999040559162, 0.999990405683675, 1, -0.418856084481767, 0.35544676217239} }; double gain=0.098531160923927; iir_t *f; f=init_iir_filter(coeff, gain, 2); return f; } /* cheby1_20 is a 3rd order low pass Chebychev1 filter with a break frequency of 5 hz and 1db of passband ripple, 20Hz sampling rate */ iir_t *init_cheby1_20() { double coeff[][6]= { {1, 1.00000959440837, 0, 1, -0.338535233784181, 0}, {1, 1.99999040559162, 0.999990405683675, 1, -0.00465798986434177, 0.60281664390126} }; double gain=0.13214070505854; iir_t *f; f=init_iir_filter(coeff, gain, 2); return f; } /* cheby1_50 is a 3rd order low pass Chebychev1 filter with a break frequency of 5 hz and 1db of passband ripple, 50Hz sampling rate */ iir_t *init_cheby1_50() { double coeff[][6]= { {1, 1.00000959440837, 0, 1, -0.723297433052849, 0}, {1, 1.99999040559162, 0.999990405683675, 1, -1.41449248783418, 0.746246746816894} }; double gain=0.0114746568820209; iir_t *f; f=init_iir_filter(coeff, gain, 2); return f; } /* The brake reaction time is defined as 70ms prior to the point at which the SVaccel dips below 0.1g */ template <class T> int findBrakingReactionFrame(timeSeries<T> SVaccel, timeSeries<T> SVspeed, timeSeries<T> throttle, timeSeries<T> brake) { const float bDelay = 0.167; // 167ms in CAMP but 70ms in virttex because // our brakes are too fast vector<int> bindx; // there can be big accel/decel spikes at the prndl transition so // start the search for the peak decel after the SV is safely in // motion bindx = find(SVspeed > 20.0/MPHPERFPS); int moving = bindx[0]; int peak = SVaccel.range(moving,SVaccel.getMaxIndex()).min(); // find the last throttle release prior to the peak decel bindx = find((throttle.range(moving,peak) > 0.05)); int beg = bindx.back(); // if this is a left-footed braker and still has the right foot on // the gas, give up and just set the window back an arbitrary 5s const float bthresh=2.0; if (brake[beg] > bthresh) { beg = SVaccel.index(SVaccel.time(beg) - 5.0); } // now find the first decel of more than 0.1 g's const double decelThresh = -0.1*32.0; bindx = find(SVaccel.range(beg,peak) < decelThresh); int j = bindx.front(); // interpolate to find the time at which SVaccel crossed thre threshold float deltaT = SVaccel.sampleTime(); float deltaA = SVaccel[j] - SVaccel[j-1]; float tcross = (deltaT/deltaA)*(decelThresh - SVaccel[j-1]) + SVaccel.time(j-1); tcross -= bDelay; return SVaccel.index(tcross); } template <class T> int findSteeringReactionFrame(timeSeries<T> SVLatAccel, timeSeries<T> swavel, timeSeries<T> swa, timeSeries<T> SVLatLanePos, timeSeries<T> SVYawRate) { const double endingLanePos = 0; const double startingLanePos = 9.84; // find the point at which the SV is 60% of the way between lanes const double Pt60 = (0.4*startingLanePos+0.6*endingLanePos); vector<int> indx; indx = find(SVLatLanePos <= Pt60); int lane60 = indx.front(); double t60 = SVLatLanePos.time(lane60); const double TswaRange=7.0; double t0 = t60 - TswaRange; // search in a time range beginning TswaRange seconds before // the 60% lane change point int minaccel = SVLatAccel.range(t0,t60).min(); double Tminaccel = SVLatAccel.time(minaccel); // smooth out the swa velocity using a 5Hz acausal filter timeSeries<T> swav_f = swavel; int sampleRate = int(floor(1.0/swavel.sampleTime()+0.5)); iir_t *lpf; if (sampleRate == 20) { lpf = init_cheby1_20(); } else { lpf = init_cheby1_50(); } filtfilt(lpf,swav_f); free_iir_filter(lpf); const double swaThresh = -6.5; //timeSeries<T> sv = swav_f.range(t0-10.0,Tminaccel); timeSeries<T> sv = swavel.range(t0-10.0,Tminaccel); timeSeries<T> sa = swa.range(t0-10.0,Tminaccel); vector<int> zerocross = find((((sv>>1)*sv)<=0.0) && (sv<=0.0) && (sa>swaThresh)); return zerocross.back(); } void outputMeasures(const string& dataFileName, const string& headerFileName, const vector<measures>& dependentMeasures, const int filenum) { // write the header file out the first time only if (filenum == 1) { ofstream headerFile(headerFileName.c_str(), ios::out); print_labels(headerFile,dependentMeasures[0]); headerFile.close(); } string measFileName = removePathExt(dataFileName) + ".dat"; ofstream mFile(measFileName.c_str(), ios::out); for (int i=0;i<dependentMeasures.size();i++) { mFile << dependentMeasures[i]; } return; } // custom state detector that finds all level changes vector<state> findAllLevels(stateFlag& sf) { vector<state> states; stateFlag::iterator<change> nextEdge; stateFlag::iterator<change> prevEdge; for (nextEdge=sf.begin<change>();nextEdge<sf.end<change>();) { prevEdge=nextEdge; state s(prevEdge.index(),(++nextEdge).index()); states.push_back(s); } return states; } // custom state detector that finds all scenarios // // a scenario is defined by a rising transition to a value n // followed by a falling edge transition to zero vector<state> findAllScenarios(stateFlag& sf) { vector<state> states; stateFlag::iterator<rising> risingEdge; stateFlag::iterator<falling> fallingEdge; stateFlag::iterator<falling> lastFallingEdge = sf.end<falling>(); for (risingEdge=sf.begin<rising>(), fallingEdge=sf.begin<falling>(); risingEdge < lastFallingEdge; ++risingEdge) { // find the first edge that falls to zero that follows // the rising edge while ((fallingEdge < risingEdge) || (*fallingEdge != 0)) { if (fallingEdge >= lastFallingEdge) { return states; } ++fallingEdge; } state s(risingEdge.index(),fallingEdge.index()-1); states.push_back(s); } return states; } string computeSubjectNo(timeSeries<float> subj) { stateFlag subState(subj); // histogram the subject numbers vector<state> subHist = findAllLevels(subState); // and then pick the subject number with the longest duration int subject; int duration = 0; for (int i=0;i<subHist.size();i++) { if ((subHist[i].length() > duration) && (subj[subHist[i].begin()] > 0)){ duration = subHist[i].length(); subject = subj[subHist[i].begin()]; } } ostringstream subjectStr; subjectStr << subject; return subjectStr.str(); } string intToCondition(const int c) { switch (c) { case 1 : return "30/30/15"; case 2 : return "30/30/28"; case 3 : return "30/30/39"; case 4 : return "45/45/15"; case 5 : return "45/45/28"; case 6 : return "45/45/39"; case 7 : return "60/60/15"; case 8 : return "60/60/28"; case 9 : return "60/60/39"; case 10 : return "30/20/0"; case 11 : return "30/10/0"; case 12 : return "60/50/0"; case 13 : return "60/30/0"; case 14 : return "60/15/0"; case 15 : return "30/0/0"; case 16 : return "45/0/0"; case 17 : return "60/0/0"; default: cout << "error...condition out of range (" << c << ")" << endl; return ""; break; } } vector<state> validScenarios(vector<state>& s, timeSeries<float>& subno, stateFlag& kin, timeSeries<int>& mode, timeSeries<int>& effort, timeSeries<int>& res, timeSeries<int>& aa) { vector<state> vs1; for (int i=0;i<s.size();i++) { bool delState = false; // if this scenario is the same as a future one // skip it for (int j=i+1;j<s.size();j++) { if (sameScenario(s[i],s[j],kin,mode,effort,res,aa)) { delState=true; break; } } // if this scenario is a test scenario // skip it if (subno[s[i].end()] <= 0) { delState=true; } if (!delState) { vs1.push_back(s[i]); } } return vs1; } bool sameScenario(state& s1, state& s2, stateFlag& kinematic_cond, timeSeries<int>& mode, timeSeries<int>& effort, timeSeries<int>& resolution, timeSeries<int>& aa) { int k1 = kinematic_cond[s1.end()]; int k2 = kinematic_cond[s2.end()]; int m1 = mode[s1.end()]; int m2 = mode[s2.end()]; int e1 = effort[s1.end()]; int e2 = effort[s2.end()]; int r1 = resolution[s1.end()]; int r2 = resolution[s2.end()]; int a1 = aa[s1.end()]; int a2 = aa[s2.end()]; bool same = (k1 == k2) && (e1 == e2) && (m1 == m2) && (r1 == r2) && (a1 == a2); return same; } float sampleTime(string s) { ifstream hfile(s.c_str()); if (!hfile.is_open()) { cout << "can't open file " << s << endl; exit(1); } float t0, tf, dt; hfile >> t0 >> tf >> dt; return dt; } string intToResolution(const int c) { string res; switch (c) { case 0: res = "1024x768"; break; case 1: res = "1600x1200"; break; case 2: res = "2048x1536"; break; default: cout << "unknown resolution: " << c << endl; res = "error"; break; } return res; } string intToAA(const int c) { string aa; switch (c) { case 0: aa = "2x2"; break; case 1: aa = "8x8"; break; default: cout << "unknown antialiasing: " << c << endl; aa = "error"; break; } return aa; }
1.6.3