24 #include "wcslib/C/wcs.h"
26 #include "wcslib/wcs.h"
36 using namespace hippodraw;
49 m_x_limits ( -180.0, +180.0 ),
50 m_y_limits ( - 90.0, + 90.0 ),
61 double xlo,
double xhi,
62 double ylo,
double yhi) :
63 BinaryTransform ( z, is_periodic, needs_grid, needs_x_ticks, needs_y_ticks ),
64 m_x_limits ( xlo, xhi ),
65 m_y_limits ( ylo, yhi ),
74 m_x_limits( t.limitX() ),
75 m_y_limits( t.limitY() ),
76 m_x_offset( t.xOffset() ),
77 m_y_offset( t.yOffset() )
148 if ( x2 < -DBL_EPSILON )
153 if ( overshoot > DBL_EPSILON ) {
163 if ( x2 < -DBL_EPSILON )
168 if ( undershoot > DBL_EPSILON)
177 if ( y2 < -DBL_EPSILON )
return moduloSubY ( y1, -y2 );
181 if ( overshoot > DBL_EPSILON ) {
192 if ( y2 < -DBL_EPSILON )
return moduloAddY ( y1, -y2 );
196 if ( undershoot > DBL_EPSILON )
return m_y_limits.
high() - undershoot;
204 double x_lo = lat.
low ();
205 double x_hi = lat.
high ();
207 double y_lo = lon.
low ();
208 double y_hi = lon.
high ();
213 double x_max = -1000;
215 double y_max = -1000;
218 double dx = ( x_hi - x_lo ) / n;
219 double dy = ( y_hi - y_lo ) / n;
223 for (
int i = 0; i <= n; i++)
228 x_min = ( x_min < x ) ? x_min : x;
229 x_max = ( x_max > x ) ? x_max : x;
230 y_min = ( y_min < y ) ? y_min : y;
231 y_max = ( y_max > y ) ? y_max : y;
235 for (
int i = 0; i <= n; i++)
240 x_min = ( x_min < x ) ? x_min : x;
241 x_max = ( x_max > x ) ? x_max : x;
242 y_min = ( y_min < y ) ? y_min : y;
243 y_max = ( y_max > y ) ? y_max : y;
247 for (
int i = 0; i <= n; i++)
252 x_min = ( x_min < x ) ? x_min : x;
253 x_max = ( x_max > x ) ? x_max : x;
254 y_min = ( y_min < y ) ? y_min : y;
255 y_max = ( y_max > y ) ? y_max : y;
259 for (
int i = 0; i <= n; i++)
264 x_min = ( x_min < x ) ? x_min : x;
265 x_max = ( x_max > x ) ? x_max : x;
266 y_min = ( y_min < y ) ? y_min : y;
267 y_max = ( y_max > y ) ? y_max : y;
270 return Rect ( x_min, y_min, x_max - x_min, y_max - y_min );
286 const vector < AxisTick > &
312 return ( (
double)abs( x - y ) <= 2.0 * ( y * FLT_EPSILON + FLT_MIN ) );
319 double rangeLength = range.
length();
322 rangeLength *= scale_factor;
325 double rmag = floor( log10( rangeLength ) );
328 double scalelow = range.
low() * scale_factor;
329 double scalehigh = range.
high() * scale_factor;
332 double pmag = max( floor( log10( abs ( scalehigh ) ) ),
333 floor( log10( abs ( scalelow ) ) ) );
351 const vector < AxisTick > &
355 double y = 0.0, ylabel;
370 bool use_pmag = abs ( pmag ) > 3.0;
388 decimals =
static_cast<int>( pmag - rmag );
391 if( !decimals ) decimals++;
408 decimals =
static_cast<int>( abs( rmag ) );
419 sprintf( pstr,
"%%1.%df", decimals );
423 double range_high = range.
high();
424 range_high *= scale_factor;
426 while( y <= range_high ||
FLT_EQUAL( range_high, y ) )
433 else if ( axistype ==
Axes::Y ) {
440 if ( use_pmag ) ylabel = value / pow( 10.0, pmag );
443 value /= scale_factor;
444 sprintf( labl, pstr, ylabel );
448 if ( tick_step == 0.0 )
break;
472 initwcs(
const std::string &transformName,
double* crpix,
473 double* crval,
double* cdelt,
double crota2,
bool galactic)
477 assert(
sizeof(wcsprm)<=2000);
483 wcsini(1, naxis,
m_wcs);
487 lon_type = (galactic?
"GLON-" :
"RA---") + transformName,
488 lat_type = (galactic?
"GLAT-" :
"DEC--") + transformName;
489 strcpy(
m_wcs->ctype[0], lon_type.c_str() );
490 strcpy(
m_wcs->ctype[1], lat_type.c_str() );
492 for(
int i=0; i<naxis; ++i){
493 m_wcs->crval[i] = crval[i];
494 m_wcs->crpix[i] = crpix[i];
495 m_wcs->cdelt[i] = cdelt[i];
507 const std::string what
508 (
"HippoDraw has not been compiled with WCSLIB support" );
509 throw std::runtime_error ( what );
523 double imgcrd[2], pixcrd[2];
524 double phi[1], theta[1];
527 double worldcrd[] ={lon,lat};
530 wcss2p(
m_wcs, ncoords, nelem, worldcrd, phi, theta, imgcrd, pixcrd, stat);
550 double worldcrd[2], imgcrd[2];
551 double phi[1], theta[1];
554 double pixcrd[] = {lon,lat};
557 wcsp2s(
m_wcs, ncoords, nelem, pixcrd, imgcrd, phi, theta, worldcrd, stat);
561 if ( returncode == 8 )
return false;
579 std::vector< double > & lat )
const
583 assert ( lat.size() == lon.size() );
584 int size = lat.size();
587 vector < double > worldcrd ( 2*size );
588 for (
int i = 0; i <
size; i++ ) {
589 worldcrd[2*i]= lon[i];
590 worldcrd[2*i+1]= lat[i];
596 vector < double > pixcrd( 2*size );
597 vector < double > imgcrd( 2*size );
598 vector < double > phi( size );
599 vector < double > theta( size );
600 vector < int > stat ( size );
603 wcss2p (
m_wcs, ncoords, nelem, &worldcrd[0], &phi[0], &theta[0],
604 &imgcrd[0], &pixcrd[0], &stat[0] );
607 for (
int i = 0; i <
size; i++ ) {
609 lat[i]=pixcrd[2*i+1];