From 67348e40e93910e7b3d1db0c3f7cefc5942b9d76 Mon Sep 17 00:00:00 2001 From: John Hill <john.r.hill@stfc.ac.uk> Date: Fri, 29 Aug 2014 11:18:33 +0100 Subject: [PATCH] Refs #10061 changes for IDF and improved loader speed --- .../inc/MantidDataHandling/Load.h | 4 +- .../inc/MantidDataHandling/LoadFITS.h | 32 +-- .../Framework/DataHandling/src/Load.cpp | 24 +- .../Framework/DataHandling/src/LoadFITS.cpp | 213 ++++++++++++------ .../inc/MantidQtMantidWidgets/MWRunFiles.ui | 3 + Code/Mantid/instrument/IMAT_Definition.xml | 60 +---- 6 files changed, 179 insertions(+), 157 deletions(-) diff --git a/Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/Load.h b/Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/Load.h index 2b0c7e6c99e..1e66e270910 100644 --- a/Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/Load.h +++ b/Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/Load.h @@ -43,8 +43,8 @@ namespace Mantid Load(); /// Algorithm's name for identification overriding a virtual method virtual const std::string name() const { return "Load"; } - ///Summary of algorithms purpose - virtual const std::string summary() const {return "Attempts to load a given file by finding an appropriate Load algorithm.";} + ///Summary of algorithms purpose + virtual const std::string summary() const {return "Attempts to load a given file by finding an appropriate Load algorithm.";} /// Algorithm's version for identification overriding a virtual method virtual int version() const { return 1; } diff --git a/Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/LoadFITS.h b/Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/LoadFITS.h index 9c54ea0a540..9d02a1a9eb6 100644 --- a/Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/LoadFITS.h +++ b/Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/LoadFITS.h @@ -79,39 +79,15 @@ namespace DataHandling /// Parses the header values for the FITS file bool parseHeader(FITSInfo &headerInfo); - void loadSingleBinFromFile(Mantid::API::MatrixWorkspace_sptr &workspace, FITSInfo &fitsInfo, MantidVecPtr &x, long spetraCount, long binIndex); + void loadChunkOfBinsFromFile(Mantid::API::MatrixWorkspace_sptr &workspace, vector<vector<double> > &yVals, vector<vector<double> > &eVals, void *&bufferAny, MantidVecPtr &x, long spetraCount, int bitsPerPixel, long binChunkStartIndex); API::MatrixWorkspace_sptr initAndPopulateHistogramWorkspace(); vector<FITSInfo> m_allHeaderInfo; - ///// Implement abstract Algorithm methods - //void init(); - ///// Implement abstract Algorithm methods - //void exec(); - - ///// Load file to a vector of strings - //void loadFile(std::string filename, std::vector<std::string>& lines); - - ///// Get Histogram type - //std::string getHistogramType(const std::vector<std::string>& lines); - - ///// Get Number of banks - //size_t getNumberOfBanks(const std::vector<std::string>& lines); - - ///// Scan imported file for bank information - //void scanBanks(const std::vector<std::string>& lines, std::vector<size_t>& bankStartIndex ); - - ///// Parse bank in file to a map - //void parseBank(std::map<std::string, double>& parammap, const std::vector<std::string>& lines, size_t bankid, size_t startlineindex, int nProf); - - ///// Find first INS line at or after lineIndex - //size_t findINSLine(const std::vector<std::string>& lines, size_t lineIndex); - - ///// Generate output workspace - //DataObjects::TableWorkspace_sptr genTableWorkspace(std::map<size_t, std::map<std::string, double> > bankparammap); - - + int m_binChunkSize; + static const int FIXED_HEADER_SIZE = 2880; + }; diff --git a/Code/Mantid/Framework/DataHandling/src/Load.cpp b/Code/Mantid/Framework/DataHandling/src/Load.cpp index ec4e7d39b89..1d701bebd5f 100644 --- a/Code/Mantid/Framework/DataHandling/src/Load.cpp +++ b/Code/Mantid/Framework/DataHandling/src/Load.cpp @@ -170,12 +170,24 @@ namespace Mantid // ... store it's name and version and check that all other files have loaders with the same name and version. std::string name = loader->name(); int version = loader->version(); + std::string ext = fileNames[0].substr(fileNames[0].find_last_of(".")); for(size_t i = 1; i < fileNames.size(); ++i) { - loader = getFileLoader(fileNames[i]); - - if( name != loader->name() || version != loader->version() ) - throw std::runtime_error("Cannot load multiple files when more than one Loader is needed."); + // If it's loading into a single file, perform a cursory check on file extensions only. + if(dynamic_cast<MultipleFileProperty*>(loader->getPointerToProperty("Filename")) != NULL) + { + if( fileNames[i].substr(fileNames[i].find_last_of(".")) != ext) + { + throw std::runtime_error("Cannot load multiple files when more than one Loader is needed."); + } + } + else + { + loader = getFileLoader(fileNames[i]); + + if( name != loader->name() || version != loader->version() ) + throw std::runtime_error("Cannot load multiple files when more than one Loader is needed."); + } } } @@ -331,7 +343,7 @@ namespace Mantid * Executes the algorithm. */ void Load::exec() - { + { std::vector<std::vector<std::string> > fileNames = getProperty("Filename"); // Test for loading as a single file @@ -351,7 +363,7 @@ namespace Mantid } void Load::loadSingleFile() - { + { std::string loaderName = getPropertyValue("LoaderName"); if( loaderName.empty() ) { diff --git a/Code/Mantid/Framework/DataHandling/src/LoadFITS.cpp b/Code/Mantid/Framework/DataHandling/src/LoadFITS.cpp index 7c2717abc86..50b9864299d 100644 --- a/Code/Mantid/Framework/DataHandling/src/LoadFITS.cpp +++ b/Code/Mantid/Framework/DataHandling/src/LoadFITS.cpp @@ -69,6 +69,7 @@ namespace DataHandling // Create FITS file information for each file selected std::vector<std::string> paths; boost::split(paths, getPropertyValue("Filename"), boost::is_any_of(",")); + m_binChunkSize = getProperty("FileChunkSize"); m_allHeaderInfo.resize(paths.size()); @@ -99,7 +100,7 @@ namespace DataHandling m_allHeaderInfo[i].timeBin = lexical_cast<double>(m_allHeaderInfo[i].headerKeys["TIMEBIN"]); m_allHeaderInfo[i].countsInImage = lexical_cast<long int>(m_allHeaderInfo[i].headerKeys["N_COUNTS"]); m_allHeaderInfo[i].numberOfTriggers = lexical_cast<long int>(m_allHeaderInfo[i].headerKeys["N_TRIGS"]); - m_allHeaderInfo[i].extension = m_allHeaderInfo[i].headerKeys["XTENSION"]; // Various extensions are available to the FITS format, and must be parsed differently if this is present + m_allHeaderInfo[i].extension = m_allHeaderInfo[i].headerKeys["XTENSION"]; // Various extensions are available to the FITS format, and must be parsed differently if this is present. Loader doesn't support this. } catch(bad_lexical_cast &) { @@ -116,6 +117,13 @@ namespace DataHandling } + // Check that the files use bit depths of either 8, 16 or 32 + if(m_allHeaderInfo[0].bitsPerPixel != 8 && m_allHeaderInfo[0].bitsPerPixel != 16 && m_allHeaderInfo[0].bitsPerPixel != 32) + { + g_log.error("FITS Loader only supports 8, 16 or 32 bits per pixel1."); + throw std::runtime_error("FITS loader only supports 8, 16 or 32 bits per pixel1."); + } + // Check the format is correct and create the Workspace if(headerValid) { @@ -137,7 +145,7 @@ namespace DataHandling // Invalid files, record error // TODO - } + } } /** @@ -156,6 +164,7 @@ namespace DataHandling // Specify as a MultipleFileProperty to alert loader we want multiple selected files to be loaded into a single workspace. declareProperty(new MultipleFileProperty("Filename", exts), "The input filename of the stored data"); + declareProperty(new PropertyWithValue<int>("FileChunkSize", 1000, Direction::Input), "Number of files to read into memory at a time - use lower values for machines with low memory"); declareProperty(new API::WorkspaceProperty<API::MatrixWorkspace>("OutputWorkspace", "", Kernel::Direction::Output)); } @@ -179,14 +188,14 @@ namespace DataHandling // Add key/values - these are separated by the = symbol. // If it doesn't have an = it's a comment to ignore. All keys should be unique - int eqPos = part.find('='); + auto eqPos = part.find('='); if(eqPos > 0) { string key = part.substr(0, eqPos); string value = part.substr(eqPos+1); // Comments are added after the value separated by a / symbol. Remove. - int slashPos = value.find('/'); + auto slashPos = value.find('/'); if(slashPos > 0) value = value.substr(0, slashPos); boost::trim(key); @@ -210,9 +219,13 @@ namespace DataHandling MantidVecPtr x; x.access().resize(m_allHeaderInfo.size() + 1); - // X = TIMEBIN value - x.access()[0] = m_allHeaderInfo[0].timeBin; - x.access()[1] = m_allHeaderInfo[0].timeBin + m_allHeaderInfo[0].tof; + // Init time bins + double binCount = 0; + for(int i=0;i<m_allHeaderInfo.size() + 1; ++i) + { + x.access()[i] = binCount; + if(i != m_allHeaderInfo.size()) binCount += m_allHeaderInfo[i].timeBin; + } long spectraCount = 0; if(m_allHeaderInfo[0].numberOfAxis > 0) spectraCount += m_allHeaderInfo[0].axisPixelLengths[0]; @@ -225,9 +238,9 @@ namespace DataHandling MatrixWorkspace_sptr retVal(new DataObjects::Workspace2D); retVal->initialize(spectraCount, m_allHeaderInfo.size()+1, m_allHeaderInfo.size()); - + IAlgorithm_sptr loadInst = createChildAlgorithm("LoadInstrument"); - + try { std::string directoryName = Kernel::ConfigService::Instance().getInstrumentDirectory(); @@ -242,10 +255,42 @@ namespace DataHandling g_log.information("Cannot load the instrument definition. " + string(ex.what()) ); } - for(int i=0; i<m_allHeaderInfo.size();++i) + int bitsPerPixel = m_allHeaderInfo[0].bitsPerPixel; // assumes all files have the same, which they should. + vector<vector<double> > yVals(spectraCount, std::vector<double>(m_binChunkSize)); + vector<vector<double> > eVals(spectraCount, std::vector<double>(m_binChunkSize)); + + // allocate memory to contain the data section of the file: + void * bufferAny = NULL; + bufferAny = malloc ((bitsPerPixel/8)*spectraCount); + + if (bufferAny == NULL) + { + fputs ("Memory error",stderr); exit (2); + } + + int steps = ceil(m_allHeaderInfo.size()/m_binChunkSize); + Progress prog(this,0.0,1.0,steps); + + + // Load a chunk of files at a time into workspace + try { - loadSingleBinFromFile(retVal, m_allHeaderInfo[i], x, spectraCount, i); + for(int i=0; i<m_allHeaderInfo.size(); i+=m_binChunkSize) + { + loadChunkOfBinsFromFile(retVal, yVals, eVals, bufferAny, x, spectraCount, bitsPerPixel, i); + prog.report(name()); + } } + catch(...) + { + // Exceptions should be handled internally, but catch here to free any memory. Belt and braces. + free(bufferAny); + g_log.error("FITS Loader unable to correctly parse files."); + throw std::runtime_error("FITS loader unable to correctly parse files."); + } + + // Memory no longer needed + free (bufferAny); retVal->mutableRun().addProperty("Filename", m_allHeaderInfo[0].filePath); @@ -270,77 +315,95 @@ namespace DataHandling - void LoadFITS::loadSingleBinFromFile(MatrixWorkspace_sptr &workspace, FITSInfo &fitsInfo, MantidVecPtr &x, long spectraCount, long binIndex) + void LoadFITS::loadChunkOfBinsFromFile(MatrixWorkspace_sptr &workspace, vector<vector<double> > &yVals, vector<vector<double> > &eVals, void *&bufferAny, MantidVecPtr &x, long spectraCount, int bitsPerPixel, long binChunkStartIndex) { - // READ DATA - ifstream istr(fitsInfo.filePath, ios::binary); - Poco::BinaryReader reader(istr); - string tmp; - reader.readRaw(2880, tmp); // Read header as full block to skip to data - // read bitdepth*naxis1 num of bits for first row? repeat naxis2 number of times. - int bytesPerRow = (fitsInfo.bitsPerPixel*fitsInfo.axisPixelLengths[0])/8; - - int allDataSizeBytes = (fitsInfo.bitsPerPixel*fitsInfo.axisPixelLengths[0]*fitsInfo.axisPixelLengths[1])/8 ; - string allData; - reader >> allData; - stringstream ss; - ss << allData; - - int currPixel, currRow = 0; - int8_t tmp8; - int16_t tmp16; - int32_t tmp32; - - for(int i=0;i<fitsInfo.axisPixelLengths[1]; ++i) // loop rows + int binsThisChunk = m_binChunkSize; + if((binChunkStartIndex + m_binChunkSize) > m_allHeaderInfo.size()) { - // Read all columns on this row. As MantidVecPt->setData expects vectors, populate data vectors for y and e. - //data.push_back(vector<double>()); - currRow = i*fitsInfo.axisPixelLengths[0]; + // No need to do extra processing if number of bins to process is lower than m_binChunkSize + binsThisChunk = m_allHeaderInfo.size() - binChunkStartIndex; + } - for(int j=0; j<fitsInfo.axisPixelLengths[0];++j) + int8_t *buffer8 = NULL; + int16_t *buffer16 = NULL; + int32_t *buffer32 = NULL; + + // create pointer of correct data type to void pointer of the buffer: + buffer8 = static_cast<int8_t*>(bufferAny); + buffer16 = static_cast<int16_t*>(bufferAny); + buffer32 = static_cast<int32_t*>(bufferAny); + + FILE *currFile = NULL; + size_t result = 0; + double val = 0; + + for(int i=binChunkStartIndex; i < binChunkStartIndex+binsThisChunk ; ++i) + { + // Read Data + currFile = fopen ( m_allHeaderInfo[i].filePath.c_str(), "rb" ); + if (currFile==NULL) { - double val = 0; - switch(fitsInfo.bitsPerPixel) - { - case 8: - ss >> tmp8; - val = static_cast<double>(tmp8); - case 16: // 2 bytes uint_16 - ss >> tmp16; - val = static_cast<double>(tmp16); - break; - case 32: - ss >> tmp32; - val = static_cast<double>(tmp32); - break; - default: - // TODO unhandled, report error. - break; - } - - currPixel = currRow + j; - workspace->setX(currPixel,x); - - //workspace->setData(currPixel,y,e); - workspace->dataY(currPixel)[binIndex] = val; - workspace->dataE(currPixel)[binIndex] = sqrt(val); - //workspace->dataX(currPixel)[binIndex] = x[0]; - - workspace->getSpectrum(currPixel)->setDetectorID(detid_t(currPixel)); - workspace->getSpectrum(currPixel)->setSpectrumNo(specid_t(currPixel+1)); + fputs ("File error",stderr); exit (1); + } + + fseek (currFile , FIXED_HEADER_SIZE , SEEK_CUR); + result = fread(bufferAny, bitsPerPixel/8, spectraCount, currFile); + + if (result != spectraCount) + { + fputs ("Reading error",stderr); exit (3); } + for(int j=0; j<spectraCount;++j) + { + if(bitsPerPixel == 8) val = static_cast<double>(buffer8[j]); + if(bitsPerPixel == 16) val = static_cast<double>(buffer16[j]); + if(bitsPerPixel == 32) val = static_cast<double>(buffer32[j]); + + yVals[j][i-binChunkStartIndex] = val; + eVals[j][i-binChunkStartIndex] = sqrt(val); + } + + // Clear memory associated with the file load + fclose (currFile); } + // Now load chunk into workspace + PARALLEL_FOR1(workspace) + for (int wi = 0; wi < spectraCount; ++wi) + { + workspace->setX(wi, x); + //MantidVecPtr y, e; + //y.access() = yVals[wi]; + //e.access() = eVals[wi]; + //retVal->setData(wi,y,e); + + //workspace.get()->dataY(wi)..insert(yVals[wi]); + + MantidVec *currY = &workspace->dataY(wi); + MantidVec *currE = &workspace->dataE(wi); + //currY->insert(currY->begin()+binChunkStartIndex, yVals[wi].begin(), yVals[wi].end()); + //currE->insert(currE->begin()+binChunkStartIndex, eVals[wi].begin(), eVals[wi].end()); + + std::copy(yVals[wi].begin(), yVals[wi].end()-(m_binChunkSize-binsThisChunk), currY->begin()+binChunkStartIndex ); + std::copy(eVals[wi].begin(), eVals[wi].end()-(m_binChunkSize-binsThisChunk), currE->begin()+binChunkStartIndex ); + //workspace->dataY(wi).insert(push_back(yVals[wi]);//.push_back(y); + //workspace->setData(wi, y, e); + // workspace->getSpectrum(currPixel)->setDetectorID(detid_t(currPixel)); + // workspace->getSpectrum(currPixel)->setSpectrumNo(specid_t(currPixel+1)); + } + } +} +} + + - } -} -} // TODO: Correctly populate X values. +// TODO: make buffer/malloc work with multiple bitdepths // about 12 seconds creating child algorithm // about 18 s loading idf @@ -349,4 +412,16 @@ namespace DataHandling // std::clock_t start; // double duration; // start = std::clock(); -//duration = ( std::clock() - start ) / (double) CLOCKS_PER_SEC; \ No newline at end of file +//duration = ( std::clock() - start ) / (double) CLOCKS_PER_SEC; + + //workspace->setData(currPixel,y,e); + //workspace->dataX(currPixel)[binIndex] = x[0]; + + + + //currPixel = currRow + j; + //workspace->dataY(currPixel)[binIndex] = val; + //workspace->dataE(currPixel)[binIndex] = sqrt(val); + + //int bytesPerRow = (m_allHeaderInfo[i].bitsPerPixel*m_allHeaderInfo[i].axisPixelLengths[0])/8; + //int allDataSizeBytes = (m_allHeaderInfo[i].bitsPerPixel*m_allHeaderInfo[i].axisPixelLengths[0]*m_allHeaderInfo[i].axisPixelLengths[1])/8 ; \ No newline at end of file diff --git a/Code/Mantid/MantidQt/MantidWidgets/inc/MantidQtMantidWidgets/MWRunFiles.ui b/Code/Mantid/MantidQt/MantidWidgets/inc/MantidQtMantidWidgets/MWRunFiles.ui index 64e18aab8e7..0e1cd0f8cea 100644 --- a/Code/Mantid/MantidQt/MantidWidgets/inc/MantidQtMantidWidgets/MWRunFiles.ui +++ b/Code/Mantid/MantidQt/MantidWidgets/inc/MantidQtMantidWidgets/MWRunFiles.ui @@ -32,6 +32,9 @@ </item> <item> <widget class="QLineEdit" name="fileEditor"> + <property name="maxLength"> + <number>500000</number> + </property> <property name="sizePolicy"> <sizepolicy hsizetype="Expanding" vsizetype="Fixed"> <horstretch>0</horstretch> diff --git a/Code/Mantid/instrument/IMAT_Definition.xml b/Code/Mantid/instrument/IMAT_Definition.xml index 6d3a86471bb..c150ddcb1a8 100644 --- a/Code/Mantid/instrument/IMAT_Definition.xml +++ b/Code/Mantid/instrument/IMAT_Definition.xml @@ -15,22 +15,11 @@ y-axis points up and the coordinate system is right handed. --> <pointing-up axis="y"/> <handedness val="right"/> </reference-frame> - <default-view axis-view="z-"/> + <default-view axis-view="z"/> </defaults> - <!-- BRIEF DESCRIPTION OF SANS2d INSTRUMENT: -Data provided by Richard Heenan (and Freddie) for the SANS2D instrument -12/06/09 this version has X & Y coords detector swapped so orientation -is correct for temporary wiring table. -18/06/09 better distances for detectors and both at L2=4m, front at X=-1.1m -26/06/09 swap front & rear as names wrong, translate front in opposite direction -21/07/09 remove the 150mm sideways shift (i.e. back to symmetrical detector coords) -to simplify manipulations in Mantid and help allow for detector mapping not quite -as expected. -01/02/10 very small chang eto pixel size 191*5.1=974.2=2*487.05 (was 487.4) -- note have to swap x= and y= in Anders output list ! -02/04/12 Put in 'no shape monitors' for possible in the future monitors -with ID 5-8 ---> + <!-- BRIEF DESCRIPTION OF IMAT INSTRUMENT: + + --> <!-- LIST OF PHYSICAL COMPONENTS (which the instrument consists of) --> <!-- source and sample-position components --> <component type="source"> @@ -47,28 +36,9 @@ with ID 5-8 </component> <type name="monitors"> <component type="monitor-tbd"> - <!-- better positions and shapes will be defined later --> <location z="7.217" name="monitor1"/> <location z="17.937" name="monitor2"/> - </component> - <component type="Moderator-Monitor3"> - <!-- transmisssion detector, either in or out of beam --> - <location z="19.497" name="monitor3"/> - </component> - <component type="monitor-tbd"> - <!-- better positions and shapes will be defined later --> - <location z="30.0" name="monitor4"/> - </component> - <!-- Putting in monitors, which are defined in raw/neuxs -files, and have detector IDs, but currently not physically present -on the instrument. Defined with no geometric shape, as they do not -physically exist, and with a dummy position --> - <component type="no shape monitor"> - <location z="0" name="placeholder monitor"/> - <location z="0" name="placeholder monitor"/> - <location z="0" name="placeholder monitor"/> - <location z="0" name="placeholder monitor"/> - </component> + </component> </type> <type name="monitor-tbd" is="monitor"> <cylinder id="some-shape"> @@ -78,22 +48,8 @@ physically exist, and with a dummy position --> <height val="0.03" /> </cylinder> </type> - <type name="Moderator-Monitor3" is="monitor"> - <percent-transparency val="99.9" /> - <cuboid id="shape"> - <left-front-bottom-point x="0.0125" y="-0.0125" z="0.0" /> - <left-front-top-point x="0.0125" y="-0.0125" z="0.005" /> - <left-back-bottom-point x="-0.0125" y="-0.0125" z="0.0" /> - <right-front-bottom-point x="0.0125" y="0.0125" z="0.0" /> - </cuboid> - <algebra val="shape" /> - </type> - <type name="no shape monitor" is="monitor" /> - <component type="detector-bank" idstart="2000000" idfillbyfirst="y" idstep="1000" idstepbyrow="1"> - <location x="1.1" z="23.281" name="front-detector"/> - </component> - <component type="detector-bank" idstart="1000000" idfillbyfirst="y" idstep="1000" idstepbyrow="1"> - <location z="23.281" name="rear-detector"/> + <component type="detector-bank" idstart="0" idfillbyfirst="y" idstep="1000" idstepbyrow="1"> + <location z="23.281" x="-0.7" name="detector"/> </component> <type name="detector-bank" is="rectangular_detector" type="pixel" xpixels="512" xstart="-0.48705" xstep="+0.0051" @@ -110,6 +66,6 @@ physically exist, and with a dummy position --> </type> <!-- DETECTOR and MONITOR ID LISTS --> <idlist idname="monitors"> - <id start="1" end="8" /> + <id start="1" end="2" /> </idlist> </instrument> \ No newline at end of file -- GitLab