More work on RenderableVolumeExpert

- Added include directory support
- Added definitions support
- Added a volume helpers
- Fixed errors in kameleon wrapper
- Added ENLIL support in kameleonwrapper
This commit is contained in:
Jonas Strandstedt
2014-04-28 15:38:40 -04:00
parent 59c92ca073
commit e4c0b88742
7 changed files with 250 additions and 37 deletions
@@ -97,6 +97,8 @@ private:
ghoul::opencl::CLWorkSize* _ws;
ghoul::filesystem::File* _kernelSourceFile;
std::vector<std::pair<ghoul::opencl::CLProgram::Option, bool> > _kernelOptions;
std::vector<std::string> _kernelIncludes;
std::vector<std::pair<std::string,std::string> > _kernelDefinitions;
bool _programUpdateOnSave;
// mutexes to prevent inconsistencies
@@ -107,6 +109,7 @@ private:
GLuint _screenQuad;
VolumeRaycasterBox* _colorBoxRenderer;
glm::vec3 _boxScaling;
};
+1
View File
@@ -49,6 +49,7 @@ public:
private:
ccmc::Model* _model;
Model _type;
ccmc::Interpolator* _interpolator;
};
+26
View File
@@ -0,0 +1,26 @@
float3 CartesianToSpherical(float3 _cartesian);
float intensityNormalizer(float intensity, float iMin, float iMax);
float3 CartesianToSpherical(float3 _cartesian) {
// Put cartesian in [-1..1] range first
_cartesian = (float3)(-1.0) + 2.0f* _cartesian;
float r = length(_cartesian);
float theta, phi;
if (r == 0.0) {
theta = phi = 0.0;
} else {
theta = acospi(_cartesian.z/r);
phi = (M_PI + atan2(_cartesian.y, _cartesian.x)) / (2.0*M_PI);
}
r = r / native_sqrt(3.0f);
// Sampler ignores w component
return (float3)(r, theta, phi);
}
float intensityNormalizer(float intensity, float iMin, float iMax) {
float i = clamp(intensity, iMin, iMax);
i = (i - iMin) / (iMax - iMin);
return clamp(i, 0.0f, 1.0f);
}
+53 -2
View File
@@ -64,7 +64,8 @@ RenderableVolumeExpert::RenderableVolumeExpert(const ghoul::Dictionary& dictiona
RenderableVolume(dictionary),
_output(nullptr),
_clBackTexture(0), _clFrontTexture(0), _clOutput(0),
_kernelSourceFile(nullptr), _programUpdateOnSave(false), _colorBoxRenderer(nullptr) {
_kernelSourceFile(nullptr), _programUpdateOnSave(false), _colorBoxRenderer(nullptr),
_boxScaling(1.0,1.0,1.0) {
_kernelLock = new std::mutex;
_textureLock = new std::mutex;
@@ -161,6 +162,30 @@ RenderableVolumeExpert::RenderableVolumeExpert(const ghoul::Dictionary& dictiona
}
}
}
ghoul::Dictionary includeDictionary;
if (kernelDictionary.hasKey("Includes") && kernelDictionary.getValue("Includes", includeDictionary)) {
auto keys = includeDictionary.keys();
for(auto key: keys) {
std::string includePath;
if(includeDictionary.getValue(key, includePath)) {
if(FileSys.directoryExists(includePath)) {
_kernelIncludes.push_back(absPath(includePath));
}
}
}
}
ghoul::Dictionary defineDictionary;
if (kernelDictionary.hasKey("Definitions") && kernelDictionary.getValue("Definitions", defineDictionary)) {
auto keys = defineDictionary.keys();
for(auto key: keys) {
std::string defintion;
if(defineDictionary.getValue(key, defintion)) {
_kernelDefinitions.push_back(std::make_pair(key, defintion));
}
}
}
}
}
if (kernelPath != "") {
@@ -168,6 +193,22 @@ RenderableVolumeExpert::RenderableVolumeExpert(const ghoul::Dictionary& dictiona
}
_colorBoxRenderer = new VolumeRaycasterBox();
double tempValue;
if(dictionary.hasKey("BoxScaling.1") && dictionary.getValue("BoxScaling.1", tempValue)) {
if(tempValue > 0.0) {
_boxScaling[0] = tempValue;
}
}
if(dictionary.hasKey("BoxScaling.2") && dictionary.getValue("BoxScaling.2", tempValue)) {
if(tempValue > 0.0) {
_boxScaling[1] = tempValue;
}
}
if(dictionary.hasKey("BoxScaling.3") && dictionary.getValue("BoxScaling.3", tempValue)) {
if(tempValue > 0.0) {
_boxScaling[2] = tempValue;
}
}
}
RenderableVolumeExpert::~RenderableVolumeExpert() {
@@ -207,6 +248,8 @@ bool RenderableVolumeExpert::initialize() {
ghoul::opengl::Texture* volume = loadVolume(_volumePaths.at(i), _volumeHints.at(i));
if(volume) {
volume->uploadTexture();
LDEBUG("Creating CL texture from GL texture with path '" << _volumePaths.at(i) << "'");
cl_mem volumeTexture = _context.createTextureFromGLTexture(CL_MEM_READ_ONLY, *volume);
_volumes.push_back(volume);
@@ -304,7 +347,7 @@ void RenderableVolumeExpert::render(const Camera *camera, const psc &thisPositio
double factor = pow(10.0,thisPosition[3]);
transform = glm::translate(transform, glm::vec3(thisPosition[0]*factor, thisPosition[1]*factor, thisPosition[2]*factor));
transform = glm::rotate(transform, time*speed, glm::vec3(0.0f, 1.0f, 0.0f));
transform = glm::scale(transform, _boxScaling);
_colorBoxRenderer->render(transform);
_textureLock->lock();
@@ -354,6 +397,14 @@ void RenderableVolumeExpert::safeKernelCompilation() {
for(auto option: _kernelOptions) {
tmpProgram.setOption(option.first, option.second);
}
for(auto defintion: _kernelDefinitions) {
tmpProgram.addDefinition(defintion.first, defintion.second);
}
// add the include directories
tmpProgram.addIncludeDirectory(_kernelIncludes);
if(tmpProgram.build()) {
ghoul::opencl::CLKernel tmpKernel = tmpProgram.createKernel("volumeraycaster");
if(tmpKernel.isValidKernel()) {
+165 -33
View File
@@ -30,12 +30,14 @@
#include <ccmc/BATSRUS.h>
#include <ccmc/ENLIL.h>
#include <iomanip>
namespace openspace {
std::string _loggerCat = "KameleonWrapper";
KameleonWrapper::KameleonWrapper(const std::string& filename, Model model) {
switch (model) {
KameleonWrapper::KameleonWrapper(const std::string& filename, Model model): _type(model) {
switch (_type) {
case Model::BATSRUS:
_model = new ccmc::BATSRUS();
if(!_model) LERROR("BATSRUS:Failed to create model instance");
@@ -53,7 +55,7 @@ KameleonWrapper::KameleonWrapper(const std::string& filename, Model model) {
if (!_interpolator) LERROR("Failed to create interpolator");
break;
default:
LERROR("Only the BATSRUS model is supported for now. Sorry.");
LERROR("No valid model type provided!");
}
}
@@ -65,47 +67,177 @@ KameleonWrapper::~KameleonWrapper() {
float* KameleonWrapper::getUniformSampledValues(const std::string& var, glm::size3_t outDimensions) {
assert(_model && _interpolator);
assert(outDimensions.x > 0 && outDimensions.y > 0 && outDimensions.z > 0);
LDEBUG("getUniformSampledValues");
assert(_type == Model::ENLIL || _type == Model::BATSRUS);
LINFO("Loading CDF data");
int size = outDimensions.x*outDimensions.y*outDimensions.z;
float* data = new float[size];
// get the grid system string
std::string gridSystem = _model->getGlobalAttribute("grid_system_1").getAttributeString();
// remove leading and trailing brackets
gridSystem = gridSystem.substr(1,gridSystem.length()-2);
// remove all whitespaces
gridSystem.erase(remove_if(gridSystem.begin(), gridSystem.end(), isspace), gridSystem.end());
// replace all comma signs with whitespaces
std::replace( gridSystem.begin(), gridSystem.end(), ',', ' ');
// tokenize
std::istringstream iss(gridSystem);
std::vector<std::string> tokens{std::istream_iterator<std::string>{iss},std::istream_iterator<std::string>{}};
// validate
if (tokens.size() != 3) {
LERROR("Something went wrong");
delete[] data;
return 0;
}
// TODO Check which coordinate system the model use. Currently assumes {x,y,z}
float xMin = _model->getVariableAttribute("x", "actual_min").getAttributeFloat();
float xMax = _model->getVariableAttribute("x", "actual_max").getAttributeFloat();
float yMin = _model->getVariableAttribute("y", "actual_min").getAttributeFloat();
float yMax = _model->getVariableAttribute("y", "actual_max").getAttributeFloat();
float zMin = _model->getVariableAttribute("z", "actual_min").getAttributeFloat();
float zMax = _model->getVariableAttribute("z", "actual_max").getAttributeFloat();
float varMin = _model->getVariableAttribute(var, "actual_min").getAttributeFloat();
float varMax = _model->getVariableAttribute(var, "actual_max").getAttributeFloat();
std::string v_x = tokens.at(0), v_y = tokens.at(1), v_z = tokens.at(2);
/*
for(auto t: tokens)
LDEBUG("t: " << t);
*/
/*
LERROR("getVariableAttributeNames");
std::vector<std::string> attributeNames = _model->getVariableAttributeNames();
for(auto name : attributeNames)
LDEBUG(name);
*/
//_model->getVa
//auto fan = std::find(attributeNames.begin(), attributeNames.end(), "");
//KameleonWrapper (Debug) grid_system_1
//KameleonWrapper (Debug) grid_1_type
LDEBUG("Using coordinate system: " << v_x << ", " << v_y << ", " << v_z);
float xMin = _model->getVariableAttribute(v_x, "actual_min").getAttributeFloat();
float xMax = _model->getVariableAttribute(v_x, "actual_max").getAttributeFloat();
float yMin = _model->getVariableAttribute(v_y, "actual_min").getAttributeFloat();
float yMax = _model->getVariableAttribute(v_y, "actual_max").getAttributeFloat();
float zMin = _model->getVariableAttribute(v_z, "actual_min").getAttributeFloat();
float zMax = _model->getVariableAttribute(v_z, "actual_max").getAttributeFloat();
float varMin = _model->getVariableAttribute(var, "actual_min").getAttributeFloat();
float varMax = _model->getVariableAttribute(var, "actual_max").getAttributeFloat();
float stepX = (xMax-xMin)/((float)outDimensions.x-1.0);
float stepY = (yMax-yMin)/((float)outDimensions.y-1.0);
float stepZ = (zMax-zMin)/((float)outDimensions.z-1.0);
// Temporary hack to keep data proportions correct, will give a lot of empty
// voxels in y and z. TODO Add spacing/voxel dimensions
float step = std::max(stepX, stepY);
step = std::max(step, stepZ);
for (int x = 0; x < outDimensions.x; ++x) {
unsigned int progress = (unsigned int)(((float)x/(float)outDimensions.x)*100.f);
if (progress % 10 == 0) {
std::cout << "Getting data from kameleon: "<< progress << "% \r" << std::flush;
float stepX = (xMax-xMin)/(static_cast<float>(outDimensions.x));
float stepY = (yMax-yMin)/(static_cast<float>(outDimensions.y));
float stepZ = (zMax-zMin)/(static_cast<float>(outDimensions.z));
LDEBUG(v_x << "Min: " << xMin);
LDEBUG(v_x << "Max: " << xMax);
LDEBUG(v_y << "Min: " << yMin);
LDEBUG(v_y << "Max: " << yMax);
LDEBUG(v_z << "Min: " << zMin);
LDEBUG(v_z << "Max: " << zMax);
LDEBUG(var << "Min: " << varMin);
LDEBUG(var << "Max: " << varMax);
int barWidth = 70;
int lastiProgress = -1;
for (int x = 0; x < outDimensions.x; ++x) {
float progress = static_cast<float>(x) / static_cast<float>(outDimensions.x-1);
int iprogress = static_cast<int>(progress*100.0f);
if (iprogress != lastiProgress) {
int pos = barWidth * progress;
int eqWidth = pos+1;
int spWidth = barWidth - pos + 2;
std::cout << "[" << std::setfill('=') << std::setw(eqWidth)
<< ">" << std::setfill(' ') << std::setw(spWidth)
<< "] " << iprogress << " % \r" << std::flush;
}
lastiProgress = iprogress;
for (int y = 0; y < outDimensions.y; ++y) {
for (int z = 0; z < outDimensions.z; ++z) {
float xPos = xMin + step*x;
float yPos = yMin + step*y;
float zPos = zMin + step*z;
int index = x + y*outDimensions.x + z*outDimensions.x*outDimensions.y;
// get interpolated data value for (xPos, yPos, zPos) and scale to [0,1]
data[index] = (_interpolator->interpolate(var, xPos, yPos, zPos)-varMin)/(varMax-varMin);
int index = x + y*outDimensions.x + z*outDimensions.x*outDimensions.y;
if(_type == Model::BATSRUS) {
float xPos = xMin + stepX*x;
float yPos = yMin + stepY*y;
float zPos = zMin + stepZ*z;
// get interpolated data value for (xPos, yPos, zPos)
float value = _interpolator->interpolate(var, xPos, yPos, zPos);
// scale to [0,1]
data[index] = (value-varMin)/(varMax-varMin);
} else if (_type == Model::ENLIL) {
//LDEBUG("data: " << theval);
// Calculate array index
//unsigned int index = r + theta*xDim_ + phi*xDim_*yDim_;
// Put r in the [0..sqrt(3)] range
float rNorm = sqrt(3.0)*(float)x/(float)(outDimensions.x-1);
// Put theta in the [0..PI] range
float thetaNorm = M_PI*(float)y/(float)(outDimensions.y-1);
// Put phi in the [0..2PI] range
float phiNorm = 2.0*M_PI*(float)z/(float)(outDimensions.z-1);
// Go to physical coordinates before sampling
float rPh = xMin + rNorm*(xMax-xMin);
float thetaPh = thetaNorm;
//phi range needs to be mapped to the slightly different
// model range to avoid gaps in the data
// Subtract a small term to avoid rounding errors when comparing
// to phiMax.
float phiPh = zMin + phiNorm/(2.0*M_PI)*(zMax-zMin-0.000001);
// Hardcoded variables (rho or rho - rho_back)
// TODO Don't hardcode, make more flexible
float varValue = 0.f;//, rho_back = 0.f, diff = 0.f;
// See if sample point is inside domain
if (rPh < xMin || rPh > xMax || thetaPh < yMin ||
thetaPh > yMax || phiPh < zMin || phiPh > zMax) {
if (phiPh > zMax) {
std::cout << "Warning: There might be a gap in the data\n";
}
// Leave values at zero if outside domain
} else { // if inside
// ENLIL CDF specific hacks!
// Convert from meters to AU for interpolator
rPh /= ccmc::constants::AU_in_meters;
// Convert from colatitude [0, pi] rad to latitude [-90, 90] degrees
thetaPh = -thetaPh*180.f/M_PI+90.f;
// Convert from [0, 2pi] rad to [0, 360] degrees
phiPh = phiPh*180.f/M_PI;
// Sample
varValue = _interpolator->interpolate(var, rPh, thetaPh, phiPh);
//rho_back = _interpolator->interpolate("rho-back",rPh,thetaPh,phiPh);
// Calculate difference (or just rho)
//diff = rho;
//diff = rho - rho_back;
// Clamp to 0
//if (diff < 0.f) diff = 0.f;
}
//if(var < 0.0f) var = 0.0f;
//data[index] = var;
data[index] = (varValue-varMin)/(varMax-varMin);
//LDEBUG("varValue:" << varValue);
//LDEBUG("data[index]:" << data[index]);
//data[index] = var;
//data[index] = diff;
}
}
}
}
std::cout << std::endl;
LINFO("Done!");
return data;
}