SyB3R - Synthetic Benchmark for 3D Reconstruction
CameraResponseCurve.h
1 #ifndef CAMERARESPONSECURVE_H
2 #define CAMERARESPONSECURVE_H
3 
4 #include <Eigen/Dense>
5 
6 #include <boost/filesystem.hpp>
7 
8 #include <vector>
9 #include <string>
10 #include <initializer_list>
11 
12 #include <iostream>
13 
14 namespace tinyxml2 {
15  class XMLDocument;
16  class XMLElement;
17 }
18 
19 
20 namespace syb3r {
21 namespace models {
22 
23 template<unsigned numLdrValues>
25 {
26  public:
27  static CameraResponseCurve<numLdrValues> readFromXML(tinyxml2::XMLElement *rootNode, const char *nodeName); // these only work for extern templates
28  void writeToXML(tinyxml2::XMLElement *rootNode, const char *nodeName) const;
29 
31  for (unsigned i = 0; i < numLdrValues; i++)
32  m_hdrValues[i] = i * 2.0f / (numLdrValues - 1);
33  }
34 
35  CameraResponseCurve(std::initializer_list<float> values) {
36  assert((values.size() == numLdrValues) && "Number of values in initializer list must match numLdrValues");
37  for (unsigned i = 0; i < numLdrValues; i++)
38  m_hdrValues[i] = values.begin()[i];
39  }
40 
41  float &ldrToHdr(unsigned ldr) { return m_hdrValues[ldr]; }
42  const float &ldrToHdr(unsigned ldr) const { return m_hdrValues[ldr]; }
43 
44  unsigned ldrFromHdr(float hdr) const;
45  protected:
46  float m_hdrValues[numLdrValues];
47 };
48 
49 extern template class CameraResponseCurve<256>;
50 
51 template<unsigned numLdrValues>
53 {
54  public:
58 
59  static CameraRGBResponseCurves<numLdrValues> readFromXML(tinyxml2::XMLElement *node); // these only work for extern templates
60  static CameraRGBResponseCurves<numLdrValues> readFromXML(const boost::filesystem::path &filename);
61  void writeToXML(tinyxml2::XMLElement *node) const;
62  void writeToXML(const boost::filesystem::path &filename) const;
63 };
64 
65 extern template class CameraRGBResponseCurves<256>;
66 
67 
68 template<unsigned numLdrValues>
69 unsigned CameraResponseCurve<numLdrValues>::ldrFromHdr(float hdr) const
70 {
71  int min = 0;
72  int max = numLdrValues-1;
73  while (min < max) {
74  int center = (max + min)/2;
75  float centerValue = m_hdrValues[center];
76  if (hdr < centerValue)
77  max = center-1;
78  else
79  min = center+1;
80  }
81  return std::max<int>(std::min<int>(min, numLdrValues-1), 0);
82 }
83 
84 template<unsigned ldrRange>
85 void computeCameraResponseFromMultipleExposures(const std::vector<unsigned> &values, const std::vector<float> &exposures, float smoothness, CameraResponseCurve<ldrRange> &curve)
86 {
87  assert(values.size() % exposures.size() == 0);
88 
89  auto computeWeight = [](unsigned v)->float{
90  if (v == 0) return 0.0f;
91  if (v == ldrRange-1) return 0.0f;
92  return 1.0;
93  };
94 
95  const unsigned numMeasurements = values.size();
96  const unsigned numPixels = values.size() / exposures.size();
97  const unsigned numCurvatures = ldrRange-2;
98  Eigen::MatrixXd designMatrix(numMeasurements + numCurvatures + 1, numPixels + ldrRange);
99  designMatrix.setZero();
100  Eigen::VectorXd designVector(numMeasurements + numCurvatures + 1);
101  designVector.setZero();
102  for (unsigned i = 0; i < numMeasurements; i++) {
103  unsigned pixel = i % numPixels;
104  unsigned image = i / numPixels;
105 
106  float weight = computeWeight(values[image * numPixels + pixel]);
107  designVector[i] = std::log(exposures[image]) * weight;
108 
109  designMatrix(i, pixel) = -1.0f * weight;
110  designMatrix(i, numPixels + values[image * numPixels + pixel]) = 1.0f * weight;
111  }
112 
113  for (unsigned i = 0; i < numCurvatures; i++) {
114  float weight = computeWeight(i+1);
115  designMatrix(numMeasurements+i, numPixels + i) = 1.0f * smoothness * weight;
116  designMatrix(numMeasurements+i, numPixels + i+1) = -2.0f * smoothness * weight;
117  designMatrix(numMeasurements+i, numPixels + i+2) = 1.0f * smoothness * weight;
118  designVector[numMeasurements+i] = 0.0f * smoothness * weight;
119  }
120 
121  designMatrix(numMeasurements+numCurvatures, numPixels + ldrRange/2) = 1.0f;
122  designVector[numMeasurements+numCurvatures] = std::log(1.0f);
123 
124  Eigen::VectorXd solution = designMatrix.jacobiSvd(Eigen::ComputeThinU | Eigen::ComputeThinV).solve(designVector);
125 
126  for (unsigned i = 0; i < ldrRange; i++)
127  curve.ldrToHdr(i) = std::exp(solution[numPixels+i]);
128 }
129 
130 #ifdef SYB3R_HAS_EXIV2
131 
132 void computeCameraResponseFromMultipleExposures(const std::vector<std::string> &jpegFilenames, unsigned numPixels, float smoothness,
134 
135 #endif
136 
138 
139 extern syb3r::models::CameraRGBResponseCurves<256> curves_Eos400D_neutral_5200K_daylight;
140 extern syb3r::models::CameraRGBResponseCurves<256> curves_Eos400D_neutral_3200K_artificial;
141 extern syb3r::models::CameraRGBResponseCurves<256> curves_Eos400D_neutral_7000K_shade;
142 
143 namespace Eos400D_neutral_3200K_artificial {
147 }
148 
149 namespace Eos400D_neutral_7000K_shade {
153 }
154 
155 
156 }
157 }
158 
159 
160 #endif // CAMERARESPONSECURVE_H
Definition: CameraPathEvaluation.cpp:10
Definition: CameraResponseCurve.h:52
Definition: CameraResponseCurve.h:24
Definition: Camera.h:6