2
$\begingroup$

NB please : executable use cases are available at the end of this question. I begin this question by showing you the problems of my program, then I explain how the latter works, and finally I end up with these both executable use cases.

What I want to do

I want to draw the bottom face of a cube. To do that, I interpolate between the 4 vertices of this face. I use two kinds of interpolation :

  1. Linear

  2. Cosine

I know it's typically the aim of GPUs, but I wanted to implement it myself to learn and try the "deep" mechanisms of 3D rendering.

Results

Each point of the following both pictures are interpolated between the 4 vertices of the bottom face (except the 8 vertices of my cube).

With linear interpolation

It's OK, no problem.

enter image description here

With cosine interpolation

Error. As you can see, there is some weird problem : face's points are inclined towards left.

enter image description here

How do I interpolate between 4 points ?

I test all the combinations of interpolation's weights and, if the sum of the weights equals 1, I draw the resulting point. I know this method is not very good (a lot of consumption of CPU and RAM), and there are P duplicates. But this problem is (1) secondary and (2) out of the range of this StackOverflow question.

In other words : the definition of the interpolation between 4 points is : P = aA + bB + cC + dD. a + b + c + d must equals 1 to draw P, otherwise it's not drawn. (a;b;c;d) \in [0;1]^4 by the way.

A concrete example is (step = 1) :

P = 0a + 0b + 0c + 0d is NOT drawn

P = 0a + 0b + 0c + 1d is drawn

P = 0a + 0b + 1c + 0d is drawn

P = 0a + 0b + 1c + 1d is NOT drawn

etc. ("I test all the combinations"). (Last is : P = 1a + 1b + 1c + 1d, which is not drawn).

To test all these combinations, I iterates on the weights, with a step, using recursion.

Case of the linear interpolation

The above definition doesn't change : P = aA + bB + cC + dD.

Case of the cosine interpolation

The above definition does change : P = aA + bB + cC + dD becomes P = a'A + b'B + c'C + d'D, with : a' = 0.5(1 - cos(a * \pi)) and the idea is the same for the three other weights.

This new definition is not random, I could explain how I found it, but it's out of the range of this topic (in résumé : a simple remap of the cosine function).

Question

Why does it works well with the linear interpolation, and not the cosine one ? Indeed, the cosine face's points are inclined towards left.

By the way : interpolations between two points work well

Linear and cosine interpolations both work very well when used between two points, as you can see below.

Linear interpolation

It's OK, no problem.

enter image description here

Cosine interpolation

It's OK, no problem.

enter image description here

Implementation and Executable

I wrote the program in Scala. First I show you the functions, then two use cases you can execute along with the given functions.

Tests all the combinations of interpolations between n points of k coordinates

def computeAllPossibleInterpolatedPoints(step : Double, points : Seq[Seq[Double]], transform: (Double) => Double) : Seq[Seq[Double]] = {
  var returned_object : Seq[Seq[Double]] = Seq.empty[Seq[Double]]
  recursiveInterpolation(0, Seq.empty[Double])

  def recursiveInterpolation(current_weight_id : Int, building_weights : Seq[Double]) : Unit = {

    (.0 to 1.0 by step).foreach(current_step => {
      if (current_weight_id < points.size - 1) {
        recursiveInterpolation(current_weight_id + 1, building_weights :+ current_step)
      } else {
        val found_solution = (building_weights :+ current_step).map(transform)
        if(BigDecimal(found_solution.sum).setScale(5, BigDecimal.RoundingMode.HALF_UP).toDouble == 1.0) {
          returned_object = returned_object :+ interpolation(found_solution, points)
        }
      }
    })
  }

  returned_object
}

Interpolates between n points of k coordinates

def interpolation(weights: Seq[Double], points: Seq[Seq[Double]], transform: (Double) => Double = null) : Seq[Double] = {
  var transformed_weights = weights
  if(transform != null) {
    transformed_weights = weights.map(transform)
  }
  if(BigDecimal(transformed_weights.sum).setScale(5, BigDecimal.RoundingMode.HALF_UP).toDouble != 1.0) {
    println("ERROR : `SUM(weights) != 1`. Returning `null`.")
    return null
  }
  if(transformed_weights.exists(weight => BigDecimal(weight).setScale(5, BigDecimal.RoundingMode.HALF_UP).toDouble < 0)
    ||
    transformed_weights.exists(weight => BigDecimal(weight).setScale(5, BigDecimal.RoundingMode.HALF_UP).toDouble > 1)) {
    println("ERROR : `EXISTS(weight) / weight < -1 OR weight > 1`. Returning `null`.")
    return null
  }

  transformed_weights.zip(points).map(
    weight_point => weight_point._2.map(coordinate => weight_point._1 * coordinate)
  ).reduce((point_a : Seq[Double], point_b : Seq[Double]) => point_a.zip(point_b).map(coordinate_points => coordinate_points._1 + coordinate_points._2))
}

Interpolation functions

def linear(weight : Double) : Double = {
  weight
}

def cosine(weight : Double) : Double = {
  (1 - Math.cos(weight * Math.PI)) * 0.5
}

Two use cases : between 2 points, each made of 2 coordinates (first use case : linear interpolation, second use case : cosine interpolation)

Note that you can add 2 new points to each of the following functions, to interpolate a face instead of a segment (for now, the below code indeed interpolates points of a segment).

val interpolated_points : Seq[Seq[Int]] = interpolator.computeAllPossibleInterpolatedPoints(0.05, Seq(Seq(22, 22), Seq(33, 33)), interpolator.linear).map(_.map(_.intValue()))

val interpolated_points : Seq[Seq[Int]] = interpolator.computeAllPossibleInterpolatedPoints(0.05, Seq(Seq(22, 22), Seq(33, 33)), interpolator.cosine).map(_.map(_.intValue()))
$\endgroup$
  • $\begingroup$ Not sure if this suits this cite (and not stackoverflow), but I definitely love the pictures. $\endgroup$ – lisyarus Jun 8 '18 at 14:30
  • $\begingroup$ Ahah thank you:)! $\endgroup$ – JarsOfJam-Scheduler Jun 8 '18 at 14:39
  • 1
    $\begingroup$ I think it can be an aliasing problem. It happens when sampling periodic functions. $\endgroup$ – Cesareo Jun 8 '18 at 14:43
  • $\begingroup$ Maybe this aliasing would occurs for linear interpolation between 2 and 4 points, and for cosine interpolation for 2 points. However, these 3 cases work well :/ In other words : why aliasing would only concern cosine interpolation between 4 points ? $\endgroup$ – JarsOfJam-Scheduler Jun 8 '18 at 15:17
  • $\begingroup$ I haven't read your implementation yet, but why do you draw $P = 0a + 0b + 0c + 0d$? According to your explanation, you should only draw points if the sum of weights is 1, but this point has weight sum 0. $\endgroup$ – N.Bach Jun 8 '18 at 15:42
3
$\begingroup$

I only remember after struggling to read through code, but I'd recommend posting pseudo-code rather than code here. Regardless, I tried to understand your implementation, and as far as I can tell, you test weights $a,b,c,d$ such that they sum to 1, and for these $a+b+c+d=1$ you directly plug them into your cosine interpolation to compute $a',b',c',d'$. In general, this will not give you a point in the face because you actually want the constraint $a'+b'+c'+d'=1$ to perform an interpolation.

Now, why does this work with a segment interpolation and not a face interpolation? The tl;dr is that $a+b=1$ will ensure $a'+b'=1$. Indeed $$\cos(a\pi)=\cos(\pi-b\pi)=-\cos(-b\pi)=-\cos(b\pi)$$ and this yields $$a'+b' = \frac 12(1-\cos(a\pi)) +\frac 12(1-\cos(b\pi))=1$$ When you have more coefficients however, this is unlikely to hold, so your interpolation is not really an interpolation.


To solve this, you can use the "transform" (linear or cosine functions) before testing the validity of the weights.

$\endgroup$
  • 1
    $\begingroup$ @JarsOfJam-Scheduler Let me rephrase my remark, it would help if your code included more comments about what you are computing. Your natural language explanation was very helpful, and it is relevant to post the actual code, but on average I'd say people here are more interested in the concepts and maths, than code. So detailed pseudo-code explaining what is computed, and in what order, is simpler to understand. Once we have figured out there are no conceptual mistakes, it is then an implementation problem. Some people prefer to help on concepts, other on implementation, so it is good to be ... $\endgroup$ – N.Bach Jun 8 '18 at 16:51
  • 1
    $\begingroup$ ...able to figure out which is which quickly. Now regarding how to correct your implementation, I suggested to use transform earlier in your code, but that implies you have to change more than just putting transform before testing your weights. If I was doing the implementation, I'd store the transformed weights, check they are valid, and then perform a "plain" interpolation using these valid weights. I don't understand the last couple lines of interpolation, but I'd get rid of map(transform) or change it to map(linear). In computeAllPossibleInterpolatedPoints I'd put ... $\endgroup$ – N.Bach Jun 8 '18 at 16:57
  • 1
    $\begingroup$ ... building_weights :+ transform(current_step) each time you sample a new weight. $\endgroup$ – N.Bach Jun 8 '18 at 16:57
  • 1
    $\begingroup$ @JarsOfJam-Scheduler That is indeed weird. It would help if you updated your code, because I cannot know if what you actually implemented is what I have in mind, or something else. Also after thinking a bit, what do you want to achieve with your cosine interpolation. Do you just want to obtain a particular distribution on each coordinates? Because in that case, this sampling method will probably not work. $\endgroup$ – N.Bach Jun 8 '18 at 19:19
  • 1
    $\begingroup$ @JarsOfJam-Scheduler I haven't fully read your update, but I tested an implementation of my own in another language. I obtained mostly the same result, and this is where I wish I had more programming experience. From what I observed, the main issue stems from testing the sum of indices being equal to 1.0. As a rule of thumb, you don't want to test equality with floating point arithmetic, especially when sampling something. Due to rounding errors in floating point arithmetic, two floats/double are rarely really equal to one another. It depends on the implementation, but ... $\endgroup$ – N.Bach Jun 9 '18 at 1:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.