programing

C: 간격으로 플로트를 감는 방법 [-pi, pi)

megabox 2023. 6. 18. 12:29
반응형

C: 간격으로 플로트를 감는 방법 [-pi, pi)

효과적으로 수행할 수 있는 멋진 C 코드를 찾고 있습니다.

while (deltaPhase >= M_PI) deltaPhase -= M_TWOPI;
while (deltaPhase < -M_PI) deltaPhase += M_TWOPI;

제가 선택할 수 있는 방법이 무엇입니까?

2013년 4월 19일 편집:

모듈로 함수가 ka.nice 및 ar_sea로 표시된 바와 같이 경계 사례를 처리하도록 업데이트되었습니다.

static const double     _PI= 3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348;
static const double _TWO_PI= 6.2831853071795864769252867665590057683943387987502116419498891846156328125724179972560696;

// Floating-point modulo
// The result (the remainder) has same sign as the divisor.
// Similar to matlab's mod(); Not similar to fmod() -   Mod(-3,4)= 1   fmod(-3,4)= -3
template<typename T>
T Mod(T x, T y)
{
    static_assert(!std::numeric_limits<T>::is_exact , "Mod: floating-point type expected");

    if (0. == y)
        return x;

    double m= x - y * floor(x/y);

    // handle boundary cases resulted from floating-point cut off:

    if (y > 0)              // modulo range: [0..y)
    {
        if (m>=y)           // Mod(-1e-16             , 360.    ): m= 360.
            return 0;

        if (m<0 )
        {
            if (y+m == y)
                return 0  ; // just in case...
            else
                return y+m; // Mod(106.81415022205296 , _TWO_PI ): m= -1.421e-14 
        }
    }
    else                    // modulo range: (y..0]
    {
        if (m<=y)           // Mod(1e-16              , -360.   ): m= -360.
            return 0;

        if (m>0 )
        {
            if (y+m == y)
                return 0  ; // just in case...
            else
                return y+m; // Mod(-106.81415022205296, -_TWO_PI): m= 1.421e-14 
        }
    }

    return m;
}

// wrap [rad] angle to [-PI..PI)
inline double WrapPosNegPI(double fAng)
{
    return Mod(fAng + _PI, _TWO_PI) - _PI;
}

// wrap [rad] angle to [0..TWO_PI)
inline double WrapTwoPI(double fAng)
{
    return Mod(fAng, _TWO_PI);
}

// wrap [deg] angle to [-180..180)
inline double WrapPosNeg180(double fAng)
{
    return Mod(fAng + 180., 360.) - 180.;
}

// wrap [deg] angle to [0..360)
inline double Wrap360(double fAng)
{
    return Mod(fAng ,360.);
}

단일 라이너 상시 솔루션:

좋아요, 두 번째 함수를 세어보면 두 개의 라이너입니다.[min,max)형태, 하지만 충분히 가깝습니다. 어쨌든 함께 병합할 수 있습니다.

/* change to `float/fmodf` or `long double/fmodl` or `int/%` as appropriate */

/* wrap x -> [0,max) */
double wrapMax(double x, double max)
{
    /* integer math: `(max + x % max) % max` */
    return fmod(max + fmod(x, max), max);
}
/* wrap x -> [min,max) */
double wrapMinMax(double x, double min, double max)
{
    return min + wrapMax(x - min, max - min);
}

그러면 간단히 사용할 수 있습니다.deltaPhase = wrapMinMax(deltaPhase, -M_PI, +M_PI).

시간입니다. 즉, 이 솔션은 이.[-PI,+PI)좋든 나쁘든

확인:

자, 저는 여러분이 제 말을 믿으리라고 기대하지 않습니다. 그래서 여기 경계 조건을 포함한 몇 가지 예가 있습니다.나는 명확성을 위해 정수를 사용하고 있지만, 그것은 거의 비슷하게 작동합니다.fmod()부동 소수점:

  • 인 것x:
    • wrapMax(3, 5) == 3:(5 + 3 % 5) % 5 == (5 + 3) % 5 == 8 % 5 == 3
    • wrapMax(6, 5) == 1:(5 + 6 % 5) % 5 == (5 + 1) % 5 == 6 % 5 == 1
  • 입니다.x:
    • 참고:이들은 정수 모듈로가 왼쪽 부호를 복사한다고 가정합니다. 그렇지 않으면 위의 ("양성") 대소문자를 받습니다.
    • wrapMax(-3, 5) == 2:(5 + (-3) % 5) % 5 == (5 - 3) % 5 == 2 % 5 == 2
    • wrapMax(-6, 5) == 4:(5 + (-6) % 5) % 5 == (5 - 1) % 5 == 4 % 5 == 4
  • 경계:
    • wrapMax(0, 5) == 0:(5 + 0 % 5) % 5 == (5 + 0) % 5 == 5 % 5 == 0
    • wrapMax(5, 5) == 0:(5 + 5 % 5) % 5 == (5 + 0) % 5== 5 % 5 == 0
    • wrapMax(-5, 5) == 0:(5 + (-5) % 5) % 5 == (5 + 0) % 5 == 5 % 5 == 0
      • 참고: 아도마.-0+0부동 소수점용.

wrapMinMax합니다: 랩핑 기은거동의일게작하다동니합: 래핑능▁function.x[min,max)하는 것과 .x - min[0,max-min) 다음에 고, (재)고, (재)고, (재)고,min결과적으로

부정적인 최대치가 어떻게 될지는 모르겠지만, 그것을 직접 확인해 보세요!

입력 각도가 임의로 높은 값에 도달할 수 있고 연속성이 중요한 경우에도 시도할 수 있습니다.

atan2(sin(x),cos(x))

이렇게 하면 특히 단일 정밀도(부동)에서 x의 높은 값에 대해 sin(x)과 cos(x)의 연속성이 모듈로보다 더 잘 보존됩니다.

실제로, exact_value_of_pi - double_value_근사 ~= 1.22e-16

반면, 대부분의 라이브러리/하드웨어는 삼각 함수를 평가할 때 모듈로를 적용하기 위해 높은 정밀도의 PI 근사치를 사용합니다(x86 계열은 다소 낮은 값을 사용하는 것으로 알려져 있음).

결과는 [-pi,pi]로 표시될 수 있습니다. 정확한 범위를 확인해야 합니다.

개인적으로, 저는 체계적으로 포장함으로써 여러 혁명에 도달하는 각도를 방지하고 부스트와 같은 fmod 솔루션을 고수할 것입니다.

▁is도 있습니다.fmod에서 합니다.math.h그러나 이 기호는 문제를 일으켜 결과를 적절한 범위에서 딱딱하게 만들기 위해 후속 작업이 필요합니다(이미 수행한 작업과 동일).다음과 같은 큰 가치의 경우deltaPhase이것은 아마도 'M_TWOPI'를 수백 번 감산/추가하는 것보다 빠를 것입니다.

deltaPhase = fmod(deltaPhase, M_TWOPI);

편집: 집중적으로 시도하지는 않았지만, 사용할 수 있을 것 같습니다.fmod이러한 방식으로 양의 값과 음의 값을 다르게 처리합니다.

    if (deltaPhase>0)
        deltaPhase = fmod(deltaPhase+M_PI, 2.0*M_PI)-M_PI;
    else
        deltaPhase = fmod(deltaPhase-M_PI, 2.0*M_PI)+M_PI;

계산 시간은 일정합니다(DeltaPhase의 절대값이 증가함에 따라 느려지는 반면 솔루션과 달리).

다음과 같은 작업을 수행할 수 있습니다.

double wrap(double x) {
    return x-2*M_PI*floor(x/(2*M_PI)+0.5);  
}

상당한 숫자 오류가 있을 것입니다.수치 오류에 대한 가장 좋은 해결책은 단계를 1/PI 또는 1/(2*PI) 단위로 확장하여 고정점으로 저장하는 것입니다.

라디안으로 작업하는 대신 1/(2θ) 축척된 각도를 사용하고 modf, floor 등을 사용합니다.라이브러리 함수를 사용하려면 라디안으로 다시 변환합니다.

이것은 또한 1만 반 회전을 회전하는 것과 1만 회전을 반 회전하는 것과 같은 효과가 있습니다. 이는 각도가 라디안일 경우에는 보장되지 않습니다. 왜냐하면 대략적인 표현을 합하는 것이 아니라 부동소수점 값에 정확한 표현을 사용하기 때문입니다.

#include <iostream>
#include <cmath>

float wrap_rads ( float r )
{
    while ( r > M_PI ) {
        r -= 2 * M_PI;
    }

    while ( r <= -M_PI ) {
        r += 2 * M_PI;
    }

    return r;
}
float wrap_grads ( float r )
{
    float i;
    r = modff ( r, &i );

    if ( r > 0.5 ) r -= 1;
    if ( r <= -0.5 ) r += 1;

    return r;
}

int main ()
{
    for (int rotations = 1; rotations < 100000; rotations *= 10 ) {
    {
        float pi = ( float ) M_PI;
        float two_pi = 2 * pi;

        float a = pi;
        a += rotations * two_pi;

        std::cout << rotations << " and a half rotations in radians " << a << " => " << wrap_rads ( a ) / two_pi << '\n' ;
    }
    {
        float pi = ( float ) 0.5;
        float two_pi = 2 * pi;

        float a = pi;
        a += rotations * two_pi;

        std::cout << rotations << " and a half rotations in grads " << a << " => " << wrap_grads ( a ) / two_pi << '\n' ;
    }
    std::cout << '\n';
}}

다음은 C++을 Boost와 함께 사용할 수 있는 이 질문을 찾는 다른 사용자를 위한 버전입니다.

#include <boost/math/constants/constants.hpp>
#include <boost/math/special_functions/sign.hpp>

template<typename T>
inline T normalizeRadiansPiToMinusPi(T rad)
{
  // copy the sign of the value in radians to the value of pi
  T signedPI = boost::math::copysign(boost::math::constants::pi<T>(),rad);
  // set the value of rad to the appropriate signed value between pi and -pi
  rad = fmod(rad+signedPI,(2*boost::math::constants::pi<T>())) - signedPI;

  return rad;
} 

C++11 버전, 부스트 종속성 없음:

#include <cmath>

// Bring the 'difference' between two angles into [-pi; pi].
template <typename T>
T normalizeRadiansPiToMinusPi(T rad) {
  // Copy the sign of the value in radians to the value of pi.
  T signed_pi = std::copysign(M_PI,rad);
  // Set the value of difference to the appropriate signed value between pi and -pi.
  rad = std::fmod(rad + signed_pi,(2 * M_PI)) - signed_pi;
  return rad;
}

두 임의의 숫자 사이에 부동 소수점 값(또는 이중)을 닫는 방법을 검색할 때 이 질문이 발생했습니다.제 경우에 대해 구체적으로 답변하지 않았기 때문에, 저는 여기에서 볼 수 있는 저만의 해결책을 생각해 보았습니다.이 값은 주어진 값을 사용하여 highperBound가 lowerBound와 완벽하게 일치하는 lowerBound와 upperBound 사이에서 동일하게 래핑합니다(즉, 360도 == 0도이므로 360도가 0으로 래핑됨).

이 대답이 좀 더 일반적인 경계 솔루션을 찾는 다른 사람들에게 도움이 되기를 바랍니다.

double boundBetween(double val, double lowerBound, double upperBound){
   if(lowerBound > upperBound){std::swap(lowerBound, upperBound);}
   val-=lowerBound; //adjust to 0
   double rangeSize = upperBound - lowerBound;
   if(rangeSize == 0){return upperBound;} //avoid dividing by 0
   return val - (rangeSize * std::floor(val/rangeSize)) + lowerBound;
}

정수에 대한 관련 질문은 여기에서 사용할 수 있습니다: C++에서 정수를 래핑하기 위한 깨끗하고 효율적인 알고리즘

[-θ, θ)에 대한 임의의 각도를 정규화하기 위해 반복적이지 않고 테스트된 두 개의 선형 솔루션:

double normalizeAngle(double angle)
{
    double a = fmod(angle + M_PI, 2 * M_PI);
    return a >= 0 ? (a - M_PI) : (a + M_PI);
}

마찬가지로 [0, 2θ)의 경우:

double normalizeAngle(double angle)
{
    double a = fmod(angle, 2 * M_PI);
    return a >= 0 ? a : (a + 2 * M_PI);
}

fmod()가 잘린 나눗셈을 통해 구현되고 배당과 같은 부호를 갖는 경우, 일반적인 문제를 다음과 같이 해결할 수 있습니다.

(-PI, PI]의 경우:

if (x > 0) x = x - 2PI * ceil(x/2PI)  #Shift to the negative regime
return fmod(x - PI, 2PI) + PI

[-PI, PI)의 경우:

if (x < 0) x = x - 2PI * floor(x/2PI)  #Shift to the positive regime
return fmod(x + PI, 2PI) - PI

[참고로 이것은 가짜 코드입니다; 제 원본은 Tcl로 쓰여졌고, 저는 그것으로 모든 사람을 고문하고 싶지 않았습니다.첫 번째 케이스가 필요했기 때문에 이 문제를 해결해야 했습니다.]

deltaPhase -= floor(deltaPhase/M_TWOPI)*M_TWOPI;

당신이 제안한 방법이 가장 좋습니다.작은 편향의 경우 가장 빠릅니다.프로그램의 각도가 계속해서 적절한 범위로 편향되는 경우에는 범위를 벗어난 큰 값에 거의 도달하지 않아야 합니다.따라서 매 라운드마다 복잡한 모듈식 산술 코드의 비용을 지불하는 것은 낭비처럼 보입니다.모듈식 산술(http://embeddedgurus.com/stack-overflow/2011/02/efficient-c-tip-13-use-the-modulus-operator-with-caution/) 에 비해 비교 비용이 저렴합니다.

C99에서:

float unwindRadians( float radians )
{
   const bool radiansNeedUnwinding = radians < -M_PI || M_PI <= radians;

   if ( radiansNeedUnwinding )
   {
      if ( signbit( radians ) )
      {
         radians = -fmodf( -radians + M_PI, 2.f * M_PI ) + M_PI;
      }
      else
      {
         radians = fmodf( radians + M_PI, 2.f * M_PI ) - M_PI;
      }
   }

   return radians;
}

glibc의 libm(newlib의 구현 포함)에 대해 링크하는 경우 __ieee754_rem_pio2f() 및 __ieee754_rem_pio2() 개인 함수에 액세스할 수 있습니다.

extern __int32_t __ieee754_rem_pio2f (float,float*);

float wrapToPI(float xf){
const float p[4]={0,M_PI_2,M_PI,-M_PI_2};

    float yf[2];
    int q;
    int qmod4;

    q=__ieee754_rem_pio2f(xf,yf);

/* xf = q * M_PI_2 + yf[0] + yf[1]                 /
 * yf[1] << y[0], not sure if it could be ignored */

    qmod4= q % 4;

    if (qmod4==2) 
      /* (yf[0] > 0) defines interval (-pi,pi]*/
      return ( (yf[0] > 0) ?  -p[2] : p[2] ) + yf[0] + yf[1];
    else
      return p[qmod4] + yf[0] + yf[1];
}

편집: 방금 당신이 libm.a에 링크해야 한다는 것을 깨달았는데, libm으로 선언된 기호를 찾을 수 없었습니다.그렇게

사용한 적이 있습니다(파이썬에서):

def WrapAngle(Wrapped, UnWrapped ):
    TWOPI = math.pi * 2
    TWOPIINV = 1.0 / TWOPI
    return  UnWrapped + round((Wrapped - UnWrapped) * TWOPIINV) * TWOPI

c-code 등가:

#define TWOPI 6.28318531

double WrapAngle(const double dWrapped, const double dUnWrapped )
{   
    const double TWOPIINV = 1.0/ TWOPI;
    return  dUnWrapped + round((dWrapped - dUnWrapped) * TWOPIINV) * TWOPI;
}

이렇게 하면 랩 도메인 +/- 2pi가 되므로 +/- pi 도메인의 경우 다음과 같이 처리해야 합니다.

if( angle > pi):
    angle -= 2*math.pi

언급URL : https://stackoverflow.com/questions/4633177/c-how-to-wrap-a-float-to-the-interval-pi-pi

반응형