Commit b2540417 authored by Zhang Xin's avatar Zhang Xin
Browse files

Initial commit

parents
#include "ephcom_wrapper.hpp"
Ephcom::Ephcom(std::string& DE_filename){
this->infp = NULL;
this->ORIG = EPHCOM_EARTH;
this->init(DE_filename);
}
Ephcom::~Ephcom(){
if( infp != NULL ){
fclose(infp);
infp = NULL;
}
}
void Ephcom::init(std::string& DE_filename){
infp = NULL;
infp = fopen(DE_filename.c_str(),"r");
if( infp == NULL ){
std::cout << __FILE__ << ":"
<< __func__ << ":L:"
<< __LINE__ << ":"
<< " failed to open: "
<< DE_filename
<< std::endl;
} else {
std::cout << "Successfully loaded ephcom data from: " << DE_filename << std::endl;
}
ephcom_readbinary_header(infp, &(this->header));
this->coords.km = 1; // not AU, use kilometers
this->coords.seconds = 0; // Timescale is days, not seconds
this->coords.bary = 1; // Center is Solar System Barycenter
this->coords.coordtype = 0; // No correction for light travel time or relativistic effects from Sun.
}
/*
targ: one of the object, SUN, ..., MERCURY (see ephcom_wrapper.hpp)
jd : Julian date
pos : position of the object
vel : velocity of the object
*/
void Ephcom::getPosVel(int targ, double jd, double pos[3], double vel[3]){
if( (targ < EPHCOM_MERCURY) || (targ > EPHCOM_SUN) ){
std::cout << __FILE__ << ":"
<< __func__ << ":L:"
<< __LINE__ << ":"
<< "targ must >= EPHCOM_MERCURY and <= EPHCOM_SUN"
<< std::endl;
}
double res[6];
this->coords.et2[0] = jd;
this->coords.et2[1] = 0.0;
ephcom_readbinary_header(this->infp, &this->header);
double *datablock = (double *) malloc(this->header.ncoeff * sizeof(double));
if ( ephcom_get_coords(this->infp, &this->header, &this->coords, datablock) == 0 ) {
ephcom_pleph(&this->coords, targ, this->ORIG, res);
pos[0] = res[0];
pos[1] = res[1];
pos[2] = res[2];
vel[0] = res[3] / 86400; // convert to km/sec
vel[1] = res[4] / 86400;
vel[2] = res[5] / 86400;
} else {
std::cout << __FILE__ << ":"
<< __func__ << ":L:"
<< __LINE__ << ":"
<< "Coordinates not found for Julian date:"
<< jd
<< std::endl;
}
free(datablock);
}
#include <iostream>
#include <iomanip>
#include <string>
#include <cstdlib>
#include "ephcom.h"
/*
#define EPHCOM_MERCURY 1
#define EPHCOM_VENUS 2
#define EPHCOM_EARTH 3
#define EPHCOM_MARS 4
#define EPHCOM_JUPITER 5
#define EPHCOM_SATURN 6
#define EPHCOM_URANUS 7
#define EPHCOM_NEPTUNE 8
#define EPHCOM_PLUTO 9
#define EPHCOM_MOON 10 // Moon, relative to Solar System center
#define EPHCOM_SUN 11
*/
#ifndef _EPHCOM_WRAPPER_HPP_
#define _EPHCOM_WRAPPER_HPP_
#define SUN EPHCOM_SUN // 太阳
#define MOON EPHCOM_MOON // 月球
#define PLUTO EPHCOM_PLUTO // 冥王星
#define NEPTUNE EPHCOM_NEPTUNE // 海王星
#define URANUS EPHCOM_URANUS // 天王星
#define SATURN EPHCOM_SATURN // 土星
#define JUPITER EPHCOM_JUPITER // 木星
#define MARS EPHCOM_MARS // 火星
#define EARTH EPHCOM_EARTH // 地球
#define VENUS EPHCOM_VENUS // 金星
#define MERCURY EPHCOM_MERCURY // 水星
class Ephcom{
private:
FILE *infp;
ephcom_Header header;
ephcom_Coords coords;
int ORIG; // origin used to calculation positions of Sun, Moon and planets. Default: EPHCOM_EARTH
public:
Ephcom(std::string& DE_filename);
~Ephcom();
void init(std::string& DE_filename);
// get position and velocity
/*
targ: one of the object, SUN, ..., MERCURY
jd : Julian date
pos : position of the object
vel : velocity of the object
*/
void getPosVel(int targ, double jd, double pos[3], double vel[3]);
};
#endif
#ifndef __EPHCOM_DLL_H
#define __EPHCOM_DLL_H
#ifdef USINGDLL
#if defined ( WIN32 )
// Visual C/C++, Borland, MinGW and Watcom
#if defined ( __VISUALC__ ) || defined ( _MSC_VER ) || defined ( __BORLANDC__ ) || defined ( __GNUC__ ) || defined ( __WATCOMC__ )
#define EPHCOMDLLEXPORT __declspec( dllexport )
#define EPHCOMDLLIMPORT __declspec( dllimport )
#else
#define EPHCOMDLLEXPORT
#define EPHCOMDLLIMPORT
#endif
#elif defined ( __CYGWIN__ )
#define EPHCOMDLLEXPORT __declspec( dllexport )
#define EPHCOMDLLIMPORT __declspec( dllimport )
#elif defined ( __GNUC__ ) && __GNUC__ > 3
// Follow ideas in http://gcc.gnu.org/wiki/Visibility for GCC version 4.x
// The following forces exported symbols specifically designated with
// EPHCOMDLLEXPORT to be visible.
#define EPHCOMDLLEXPORT __attribute__ ( ( visibility( "default" ) ) )
#define EPHCOMDLLIMPORT
#endif
#endif
// For an unknown compiler or static build we clear the macros
#ifndef EPHCOMDLLEXPORT
#define EPHCOMDLLEXPORT
#define EPHCOMDLLIMPORT
#endif
// The IMPEXP macros will always be set to DLLIMPORT (even for
// the static library, but DLLIMPORT is empty in this case), if
// cmake didn't set the corresponding macro xxxx_EXPORTS when the
// corresponding library is built (DLLIMPEXP is set to DLLEXPORT
// then)
#if defined ( gnulliver_EXPORTS )
#define GNULLIVERDLLIMPEXP EPHCOMDLLEXPORT
#define GNULLIVERDLLIMPEXP_DATA( type ) EPHCOMDLLEXPORT type
#else
#define GNULLIVERDLLIMPEXP EPHCOMDLLIMPORT
#define GNULLIVERDLLIMPEXP_DATA( type ) EPHCOMDLLIMPORT type
#endif
#if defined ( ephcom_EXPORTS )
#define EPHCOMDLLIMPEXP EPHCOMDLLEXPORT
#define EPHCOMDLLIMPEXP_DATA( type ) EPHCOMDLLEXPORT type
#else
#define EPHCOMDLLIMPEXP EPHCOMDLLIMPORT
#define EPHCOMDLLIMPEXP_DATA( type ) EPHCOMDLLIMPORT type
#endif
#if defined ( ephcomfc_EXPORTS )
#define EPHCOMFCDLLIMPEXP EPHCOMDLLEXPORT
#define EPHCOMFCDLLIMPEXP_DATA( type ) EPHCOMDLLEXPORT type
#else
#define EPHCOMFCDLLIMPEXP EPHCOMDLLIMPORT
#define EPHCOMFCDLLIMPEXP_DATA( type ) EPHCOMDLLIMPORT type
#endif
#endif // __EPHCOM_DLL_H
//! @file
//! Source code for the gnulliver library.
//!
// gnulliver.c - Travels the middle of the road to swiftly gnullify
// the Big-endian / Little-endian controversy.
// Only works on machines with IEEE floating point.
// Copyright (C) 1994-2004 Paul Hardy
// Copyright (C) 2011 Alan W. Irwin
//
// This file is part of the timeephem software project.
//
// timeephem is free software; you can redistribute it and/or modify
// it under the terms of the GNU Library General Public License as published
// by the Free Software Foundation; version 2 of the License.
//
// timeephem is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU Library General Public License for more details.
//
// You should have received a copy of the GNU Library General Public License
// along with timeephem; if not, write to the Free Software
// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
// Swaps values between host order and network order without calls
// to htonl(), etc. This is a symmetric conversion, so the same
// function works in both directions. For example, to convert an
// integer between host byte ordering and network byte ordering (in
// either direction), use
//
// int i, j;
// j = gnulliver32(i);
//
// then use j in an equation (if you just read i from the net), or
// write j to the net (if i is an integer on your machine).
//
// The gnulliver routines work by setting a double to a certain value,
// then examining a unioned byte array.
//
// The functions are:
//
// char gnulliver() - determine (or recalculate) endianness of host:
// GNULLIVER_BIG, GNULLIVER_LITTLE, or GNULLIVER_DEC, or whether to
// just do no swapping at all (GNULLIVER_NATIVE) to read and
// write native endian result.
// unsigned short gnulliver16(unsigned short input)
// unsigned gnulliver32(unsigned input)
// unsigned long gnulliver64(unsigned long input)
// unsigned long long gnulliver128(unsigned long long input)
// float gnulliver32f(float input)
// double gnulliver64f(double input)
// long double gnulliver128f(long double input)
//
// The code in the routines might appear largely redundant, and it is.
// Loops are also unrolled. All this should minimize function calls and
// maximize speed.
//
#include "gnulliver.h"
//
// GNULLIVER_TEST is the magic number we use to determine byte swapping in
// a 64-bit word. It assumes IEEE 754 floating point, and encodes to
// hexadecimal 0x40302010 00000000 on a Big-endian (network order) machine.
// The mantissa is large enough that it might work on a non-IEEE 754 compatible
// machine just by coincidence, but no guarantees.
//
// gnulliver currently assumes that if we swap 32-bit values in a 64-bit word,
// then we also swap 64-bit values in a 128-bit word (long double), which is
// true if the machine is Little-endian. If this does not hold true on a
// particular platform, let me know and I'll modify the code.
//
#define GNULLIVER_TEST ( 16.0 + ( 1.0 / 8.0 + 1.0 / 4096.0 ) )
//! Discover the endianess of the host.
//!
//! @returns GNULLIVER_LITTLE, GNULLIVER_DEC, GNULLIVER_BIG depending on
//! whether the host uses little-, DEC-, or big-endian order for its bytes.
//! This information is used to swap bytes between host byte order and
//! assumed network order (big-endian) for the binary file.
//! Alternatively, gnulliver() returns GNULLIVER_NATIVE for the case
//! where no byte swaps are desired for i/o on a binary file.
//!
unsigned char
gnulliver()
{
#ifdef NO_GNULLIVER_SWAP
return ( GNULLIVER_NATIVE );
#else
unsigned char result;
union
{
double d;
unsigned char ch[8];
} bytes;
bytes.d = GNULLIVER_TEST;
result = 0;
if ( bytes.ch[0] == 0 ) // swap 32-bit words in a 64-bit number
{
if ( bytes.ch[4] == 0x10 ) // Assume 0x00000000 10203040
result = GNULLIVER_LITTLE;
else // Assume 0x00000000 20104030
result = GNULLIVER_DEC;
}
else // Assume 0x40302010 00000000
result = GNULLIVER_BIG;
return ( result );
#endif
}
//! Optional byteswap of short depending on gnulliver() result.
//!
//! @param input [IN ONLY]Value of integer to be optionally byteswapped.
//! @returns (optionally) byteswapped integer.
//!
unsigned short
gnulliver16( unsigned short input )
{
static unsigned endian, init = 1;
unsigned char tmpch;
union
{
unsigned short u;
unsigned char ch[2];
} bytes;
if ( init )
{
endian = gnulliver();
init = 0;
}
bytes.u = input;
if ( endian != GNULLIVER_BIG && endian != GNULLIVER_NATIVE )
{
tmpch = bytes.ch[0];
bytes.ch[0] = bytes.ch[1]; bytes.ch[1] = tmpch;
}
return ( bytes.u );
}
//! Optional byteswap of integer depending on gnulliver() result.
//!
//! @param input [IN ONLY]Value of integer to be optionally byteswapped.
//! @returns (optionally) byteswapped integer.
//!
unsigned
gnulliver32( unsigned input )
{
static unsigned endian, init = 1;
char tmpch;
union
{
unsigned u;
unsigned char ch[4];
} bytes;
if ( init )
{
endian = gnulliver();
init = 0;
}
bytes.u = input;
if ( endian != GNULLIVER_BIG && endian != GNULLIVER_NATIVE )
{
if ( endian == GNULLIVER_LITTLE )
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[3]; bytes.ch[3] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[2]; bytes.ch[2] = tmpch;
}
else // endian == GNULLIVER_DEC
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[2]; bytes.ch[2] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[3]; bytes.ch[3] = tmpch;
}
}
return ( bytes.u );
}
//! Optional byteswap of long depending on gnulliver() result.
//!
//! @param input [IN ONLY]Value of integer to be optionally byteswapped.
//! @returns (optionally) byteswapped integer.
//!
unsigned long
gnulliver64( unsigned long input )
{
static unsigned endian, init = 1;
char tmpch;
union
{
unsigned long u;
char ch[8];
} bytes;
if ( init )
{
endian = gnulliver();
init = 0;
}
bytes.u = input;
if ( endian != GNULLIVER_BIG && endian != GNULLIVER_NATIVE )
{
if ( endian == GNULLIVER_LITTLE )
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[7]; bytes.ch[7] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[6]; bytes.ch[6] = tmpch;
tmpch = bytes.ch[2]; bytes.ch[2] = bytes.ch[5]; bytes.ch[5] = tmpch;
tmpch = bytes.ch[3]; bytes.ch[3] = bytes.ch[4]; bytes.ch[4] = tmpch;
}
else // Assume GNULLIVER_DEC
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[6]; bytes.ch[6] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[7]; bytes.ch[7] = tmpch;
tmpch = bytes.ch[2]; bytes.ch[2] = bytes.ch[4]; bytes.ch[4] = tmpch;
tmpch = bytes.ch[3]; bytes.ch[3] = bytes.ch[5]; bytes.ch[5] = tmpch;
}
}
return ( bytes.u );
}
//! Optional byteswap of long long depending on gnulliver() result.
//!
//! @param input [IN ONLY]Value of integer to be optionally byteswapped.
//! @returns (optionally) byteswapped integer.
//!
unsigned long long
gnulliver128( unsigned long long input )
{
static unsigned endian, init = 1;
unsigned char tmpch;
union
{
unsigned long long u;
char ch[16];
} bytes;
if ( init )
{
endian = gnulliver();
init = 0;
}
bytes.u = input;
if ( endian != GNULLIVER_BIG && endian != GNULLIVER_NATIVE )
{
if ( endian == GNULLIVER_LITTLE )
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[15]; bytes.ch[15] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[14]; bytes.ch[14] = tmpch;
tmpch = bytes.ch[2]; bytes.ch[2] = bytes.ch[13]; bytes.ch[13] = tmpch;
tmpch = bytes.ch[3]; bytes.ch[3] = bytes.ch[12]; bytes.ch[12] = tmpch;
tmpch = bytes.ch[4]; bytes.ch[4] = bytes.ch[11]; bytes.ch[11] = tmpch;
tmpch = bytes.ch[5]; bytes.ch[5] = bytes.ch[10]; bytes.ch[10] = tmpch;
tmpch = bytes.ch[6]; bytes.ch[6] = bytes.ch[ 9]; bytes.ch[ 9] = tmpch;
tmpch = bytes.ch[7]; bytes.ch[7] = bytes.ch[ 8]; bytes.ch[ 8] = tmpch;
}
else // Assume GNULLIVER_DEC
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[14]; bytes.ch[14] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[15]; bytes.ch[15] = tmpch;
tmpch = bytes.ch[2]; bytes.ch[2] = bytes.ch[12]; bytes.ch[12] = tmpch;
tmpch = bytes.ch[3]; bytes.ch[3] = bytes.ch[13]; bytes.ch[13] = tmpch;
tmpch = bytes.ch[4]; bytes.ch[4] = bytes.ch[10]; bytes.ch[10] = tmpch;
tmpch = bytes.ch[5]; bytes.ch[5] = bytes.ch[11]; bytes.ch[11] = tmpch;
tmpch = bytes.ch[6]; bytes.ch[6] = bytes.ch[ 8]; bytes.ch[ 8] = tmpch;
tmpch = bytes.ch[7]; bytes.ch[7] = bytes.ch[ 9]; bytes.ch[ 9] = tmpch;
}
}
return ( bytes.u );
}
//! Optional byteswap of float depending on gnulliver() result.
//!
//! @param input [IN ONLY]Value of floating-point variable to be
//! optionally byteswapped.
//! @returns (optionally) byteswapped floating point value.
//!
float
gnulliver32f( float input )
{
static unsigned endian, init = 1;
unsigned char tmpch;
union
{
float x;
unsigned char ch[4];
} bytes;
if ( init )
{
endian = gnulliver();
init = 0;
}
bytes.x = input;
if ( endian != GNULLIVER_BIG && endian != GNULLIVER_NATIVE )
{
if ( endian == GNULLIVER_LITTLE )
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[3]; bytes.ch[3] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[2]; bytes.ch[2] = tmpch;
}
else // Assume GNULLIVER_DEC
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[2]; bytes.ch[2] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[3]; bytes.ch[3] = tmpch;
}
}
return ( bytes.x );
}
//! Optional byteswap of double depending on gnulliver() result.
//!
//! @param input [IN ONLY]Value of floating-point variable to be
//! optionally byteswapped.
//! @returns (optionally) byteswapped floating point value.
//!
double
gnulliver64f( double input )
{
static unsigned endian, init = 1;
unsigned char tmpch;
union
{
double x;
unsigned char ch[8];
} bytes;
if ( init )
{
endian = gnulliver();
init = 0;
}
bytes.x = input;
if ( endian != GNULLIVER_BIG && endian != GNULLIVER_NATIVE )
{
if ( endian == GNULLIVER_LITTLE )
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[7]; bytes.ch[7] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[6]; bytes.ch[6] = tmpch;
tmpch = bytes.ch[2]; bytes.ch[2] = bytes.ch[5]; bytes.ch[5] = tmpch;
tmpch = bytes.ch[3]; bytes.ch[3] = bytes.ch[4]; bytes.ch[4] = tmpch;
}
else // Assume GNULLIVER_DEC
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[6]; bytes.ch[6] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[7]; bytes.ch[7] = tmpch;
tmpch = bytes.ch[2]; bytes.ch[2] = bytes.ch[4]; bytes.ch[4] = tmpch;
tmpch = bytes.ch[3]; bytes.ch[3] = bytes.ch[5]; bytes.ch[5] = tmpch;
}
}
return ( bytes.x );
}
//! Optional byteswap of long double depending on gnulliver() result.
//! Note: This function's result will be unpredictable if the type
//! long double is not a 128-bit word.
//!
//! @param input [IN ONLY]Value of floating-point variable to be
//! optionally byteswapped.
//! @returns (optionally) byteswapped floating point value.
//!
long double
gnulliver128f( long double input )
{
static unsigned endian, init = 1;
unsigned char tmpch;
union
{
long double x;
unsigned char ch[16];
} bytes;
if ( init )
{
endian = gnulliver();
init = 0;
}
bytes.x = input;
if ( endian != GNULLIVER_BIG && endian != GNULLIVER_NATIVE )
{
if ( endian == GNULLIVER_LITTLE )
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[15]; bytes.ch[15] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[14]; bytes.ch[14] = tmpch;
tmpch = bytes.ch[2]; bytes.ch[2] = bytes.ch[13]; bytes.ch[13] = tmpch;
tmpch = bytes.ch[3]; bytes.ch[3] = bytes.ch[12]; bytes.ch[12] = tmpch;
tmpch = bytes.ch[4]; bytes.ch[4] = bytes.ch[11]; bytes.ch[11] = tmpch;
tmpch = bytes.ch[5]; bytes.ch[5] = bytes.ch[10]; bytes.ch[10] = tmpch;
tmpch = bytes.ch[6]; bytes.ch[6] = bytes.ch[ 9]; bytes.ch[ 9] = tmpch;
tmpch = bytes.ch[7]; bytes.ch[7] = bytes.ch[ 8]; bytes.ch[ 8] = tmpch;
}
else // Assume GNULLIVER_DEC
{
tmpch = bytes.ch[0]; bytes.ch[0] = bytes.ch[14]; bytes.ch[14] = tmpch;
tmpch = bytes.ch[1]; bytes.ch[1] = bytes.ch[15]; bytes.ch[15] = tmpch;
tmpch = bytes.ch[2]; bytes.ch[2] = bytes.ch[12]; bytes.ch[12] = tmpch;
tmpch = bytes.ch[3]; bytes.ch[3] = bytes.ch[13]; bytes.ch[13] = tmpch;
tmpch = bytes.ch[4]; bytes.ch[4] = bytes.ch[10]; bytes.ch[10] = tmpch;
tmpch = bytes.ch[5]; bytes.ch[5] = bytes.ch[11]; bytes.ch[11] = tmpch;
tmpch = bytes.ch[6]; bytes.ch[6] = bytes.ch[ 8]; bytes.ch[ 8] = tmpch;
tmpch = bytes.ch[7]; bytes.ch[7] = bytes.ch[ 9]; bytes.ch[ 9] = tmpch;
}
}
return ( bytes.x );
}
//! Optional swap of char array elements depending on gnulliver()
//! result.
//!
//! @param ch [IN AND OUT ]Pointer to char array containing 8 elements
//! which upon return will be optionally swapped.
//! @returns pointer to char array which contains the (optionally)
//! swapped elements.
//!
unsigned char *
gnulliver64c( unsigned char *ch )
{
static unsigned endian, init = 1;
unsigned char tmpch;
if ( init )
{
endian = gnulliver();
init = 0;
}
if ( endian != GNULLIVER_BIG && endian != GNULLIVER_NATIVE )
{
if ( endian == GNULLIVER_LITTLE )
{
tmpch = ch[0]; ch[0] = ch[7]; ch[7] = tmpch;
tmpch = ch[1]; ch[1] = ch[6]; ch[6] = tmpch;
tmpch = ch[2]; ch[2] = ch[5]; ch[5] = tmpch;
tmpch = ch[3]; ch[3] = ch[4]; ch[4] = tmpch;
}
else // Assume GNULLIVER_DEC
{
tmpch = ch[0]; ch[0] = ch[6]; ch[6] = tmpch;
tmpch = ch[1]; ch[1] = ch[7]; ch[7] = tmpch;
tmpch = ch[2]; ch[2] = ch[4]; ch[4] = tmpch;
tmpch = ch[3]; ch[3] = ch[5]; ch[5] = tmpch;
}
}
return ( ch );
}
// gnulliver.h - header information for the gnulliver library.
// Copyright (C) 1994-2004 Paul Hardy
// Copyright (C) 2011 Alan W. Irwin
//
// This file is part of the timeephem software project.
//
// timeephem is free software; you can redistribute it and/or modify
// it under the terms of the GNU Library General Public License as published
// by the Free Software Foundation; version 2 of the License.
//
// timeephem is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU Library General Public License for more details.
//
// You should have received a copy of the GNU Library General Public License
// along with timeephem; if not, write to the Free Software
// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
#ifndef __GNULLIVER_H__
#define __GNULLIVER_H__
//
// Travels the middle of the road to swiftly gnullify
// the Big-endian / Little-endian controversy.
// Only works on machines with IEEE floating point.
//
// Swaps values between host order and network order without calls
// to htonl(), etc. This is a symmetric conversion, so the same
// function works in both directions. For example, to convert an
// integer between host byte ordering and network byte ordering (in
// either direction), use
//
// int i, j;
// j = gnulliver32(i);
//
// then use j in an equation (if you just read i from the net), or
// write j to the net (if i is an integer on your machine).
//
// The gnulliver routines work by setting a double to a certain value,
// then examining a unioned byte array.
//
// The functions are:
//
// char gnulliver() - determine (or recalculate) endianness of host:
// GNULLIVER_BIG, GNULLIVER_LITTLE, GNULLIVER_DEC
// unsigned short gnulliver16(unsigned short input)
// unsigned gnulliver32(unsigned input)
// unsigned long gnulliver64(unsigned long input)
// unsigned long long gnulliver128(unsigned long long input)
// float gnulliver32f(float input)
// double gnulliver64f(double input)
// long double gnulliver128f(long double input)
//
// Set up for dealing with function visibility issues.
#include "ephcomdll.h"
#define GNULLIVER_VERSION "1.0"
#define GNULLIVER_SWAP16 4 // Swap adjacent bytes in 2-byte string
#define GNULLIVER_SWAP32 8 // Swap 2-byte pairs in 4-byte string
#define GNULLIVER_SWAP64 16 // Swap 4-byte pairs in 8-byte string
#define GNULLIVER_SWAP128 32 // Swap 8-byte pairs in 16-byte string
//
// What we have to do to convert between network (Big-Endian) order and
// host machine order.
//
#define GNULLIVER_NATIVE 1 // No swapping regardless.
#define GNULLIVER_BIG 2 // No swapping if this machine is Big-endian
#define GNULLIVER_LITTLE ( GNULLIVER_SWAP16 | GNULLIVER_SWAP32 | GNULLIVER_SWAP64 | GNULLIVER_SWAP128 )
// This should work on a PDP-11 or VAX, but I have no way to test it.
#define GNULLIVER_DEC ( GNULLIVER_SWAP32 | GNULLIVER_SWAP64 | GNULLIVER_SWAP128 )
// Public function declarations for the gnulliver library.
GNULLIVERDLLIMPEXP unsigned char
gnulliver();
GNULLIVERDLLIMPEXP unsigned short
gnulliver16( unsigned short input );
GNULLIVERDLLIMPEXP unsigned
gnulliver32( unsigned input );
GNULLIVERDLLIMPEXP unsigned long
gnulliver64( unsigned long input );
GNULLIVERDLLIMPEXP unsigned char *
gnulliver64c( unsigned char *ch );
GNULLIVERDLLIMPEXP unsigned long long
gnulliver128( unsigned long long input );
GNULLIVERDLLIMPEXP float
gnulliver32f( float input );
GNULLIVERDLLIMPEXP double
gnulliver64f( double input );
GNULLIVERDLLIMPEXP long double
gnulliver128f( long double input );
#endif // __GNULLIVER_H__
#include <iostream>
#include"StrayLight.h"
using namespace std;
int main()
{
double ju=2.4608437604166665e+06;//观测时的儒略时
//以下所有位置和方向均为赤道坐标系三维坐标
double sat[3] = {5873.752, -1642.066, 2896.744};//卫星所处位置
double ob[3]={0.445256,0.76061,-0.47246};//观测方向
//double py[3]={0.484277, 0.239253,0.841565};//某合理Y轴
double py1[3],py2[3];//未知的Y轴
ComposeY(ob,py1,py2);//计算Y轴
double E[7];//7谱段的辐照度
Init();//初始化
Calculate(ju,sat,ob,py1,E);//计算
Calculate(ju,sat,ob,py2,E);
return 0;
}
# ra dec lon(ecliptic) lat(ecliptic) pos_angle time(julian) sat_x sat_y sat_z sun_x sun_y sun_z moon_x moon_y moon_z sat_vx sat_vy sat_vz exp_time isDeep
244.972773 39.895901 229.6933 59.8819 109.59452000342756 2460539.75361268 -2218.3663 -5817.8559 2597.9467 -123291049.6698 80698075.2698 34981832.8099 129217.8762 -306193.0456 -168270.0801 5.058199026407237 -3.85818799871231 -4.322908628205369 150.0 -1.0
245.19954 39.834237 230.0613 59.8819 109.42155570115811 2460539.75590434 -1169.3163 -6427.7801 1683.3766 -123294376.7625 80693691.118 34979932.4325 129409.3308 -306117.8806 -168231.5896 5.495628791263698 -2.2766122251687193 -4.877784128846088 150.0 -1.0
245.376849 40.538238 229.8433 60.5906 109.72521611662427 2460543.53323198 -4088.9785 -4865.03 2348.6771 -128524732.1114 73309196.4688 31779164.5982 349033.2322 -76382.1208 -45844.3874 3.513272716954816 -5.1280401490294025 -4.4997767748736806 150.0 -1.0
245.939689 40.198882 230.894 60.4134 109.17526390397103 2460543.53552584 -3296.1182 -5750.7175 1405.5287 -128527752.3153 73304618.787 31777180.4968 349086.0107 -76197.4206 -45744.3873 4.456135176128555 -3.773136335412346 -4.9799999694710095 150.0 -1.0
246.809804 39.784988 232.438 60.2362 108.39295915033918 2460543.53783608 -2330.5509 -6351.0086 384.2962 -128530793.8798 73300008.3241 31775182.1872 349139.047 -76011.3767 -45643.6581 5.179717934292967 -2.216473482430956 -5.210878017528557 150.0 -1.0
247.037801 39.727005 232.8099 60.2362 108.21592927146807 2460543.54012774 -1255.3999 -6627.1633 -648.2573 -128533810.807 73295434.8053 31773199.8903 349191.5393 -75826.8029 -45543.7232 5.637231322125444 -0.5605653312886716 -5.176883766496303 150.0 -1.0
247.26584 39.669564 233.1818 60.2362 108.038749009507 2460543.54241941 -117.3694 -6571.3891 -1648.2421 -128536827.5433 73290861.1758 31771217.5454 349243.9141 -75642.204 -45443.7731 5.812145421308514 1.122583740434493 -4.883329814897479 150.0 -1.0
245.254413 40.383388 229.7713 60.4134 109.70843030522704 2460543.59745489 -4101.478 -4866.9578 2323.1349 -128609218.3748 73180990.3315 31723596.3195 350466.4112 -71201.6169 -43038.9593 3.509223042501617 -5.1159759164056595 -4.516449863272101 150.0 -1.0
245.898224 40.022321 230.9504 60.2362 109.09954624314884 2460543.599753 -3307.4819 -5751.6002 1375.5044 -128612238.7884 73176401.0903 31721607.2092 350515.9823 -71015.8914 -42938.3589 4.456793269669106 -3.758275091286123 -4.990470905338498 150.0 -1.0
246.912851 39.57222 232.7235 60.0591 108.21113555960895 2460543.60207086 -2338.0901 -6350.3531 349.5009 -128615284.9693 73171772.2898 31719600.9527 350565.8595 -70828.5455 -42836.8793 5.185149246981382 -2.196511614110932 -5.213861687802819 150.0 -1.0
247.140446 39.514777 233.0935 60.0591 108.03541662103602 2460543.60436253 -1261.6936 -6622.571 -682.7541 -128618296.5268 73167195.698 31717617.3251 350615.0543 -70643.2938 -42736.5323 5.644314272476549 -0.5413230985595874 -5.171075108238483 150.0 -1.0
246.685297 39.630202 232.3535 60.0591 108.38670333146824 2460543.6066542 -122.112 -6563.1261 -1680.7261 -128621307.8934 73162618.9958 31715633.6496 350664.1312 -70458.0187 -42636.1712 5.820605356021247 1.1401236903839163 -4.869028668063038 150.0 -1.0
246.561326 39.475239 232.2699 59.8819 108.38013823250273 2460543.60894586 1023.5967 -6175.0802 -2594.3 -128624319.0682 73158042.184 31713649.9267 350713.0903 -70272.7204 -42535.7958 5.7056975349860295 2.764332291953906 -4.322172590177615 150.0 -1.0
245.623258 40.847614 229.9894 60.9449 109.75838193462751 2460543.6682487 -1437.2044 -6598.2001 -561.399 -128702174.3868 73039567.2025 31662299.4501 351938.9714 -65469.7753 -39933.5032 5.599335756251094 -0.7771821850383276 -5.189864591059177 150.0 -1.0
245.499769 40.692917 229.916 60.7677 109.74184789995181 2460543.67054037 -302.0605 -6585.6621 -1566.1018 -128705180.4091 73034987.4487 31660314.4531 351984.7548 -65283.8767 -39832.7597 5.820956241348654 0.9040687362721656 -4.91795858765272 150.0 -1.0
245.605624 40.476178 230.2194 60.5906 109.54618182037997 2460543.85286552 -4551.9642 -4177.9517 2779.7576 -128943724.7839 72670271.9443 31502236.3919 355246.6721 -50430.2867 -31778.0894 2.7806315862835618 -5.809120183557752 -4.171875803399416 150.0 -1.0
246.168179 40.138468 231.2682 60.4134 108.9972405764492 2460543.85515932 -3891.5017 -5214.6838 1890.235 -128946718.1231 72665679.1318 31500245.7379 355282.9058 -50242.6873 -31676.2957 3.858369242079789 -4.610655120031879 -4.76915223102651 150.0 -1.0
246.581847 39.843512 232.0661 60.2362 108.5698372158382 2460543.85745099 -3036.7067 -5989.2387 906.5607 -128949708.4894 72661090.477 31498256.886 355318.9862 -50055.2459 -31574.5862 4.741932037004517 -3.1813221279298887 -5.127064128400889 150.1976 -1.0
247.368081 39.457872 233.4635 60.0591 107.85954794587214 2460543.85975863 -2022.3883 -6466.0294 -129.8582 -128952719.5082 72656469.7232 31496254.1216 355355.1975 -49866.4812 -31472.1571 5.39186183170159 -1.5809240736052743 -5.2270457437822415 150.0 -1.0
247.49392 39.612666 233.5537 60.2362 107.86141979573472 2460543.8620503 -913.3254 -6614.5703 -1152.6821 -128955709.4883 72651880.8508 31494265.1755 355391.0381 -49679.0071 -31370.4268 5.766144793522017 0.08734229075344047 -5.06308085288083 150.0 -1.0
247.540067 39.78741 233.5135 60.4134 107.92582448507959 2460543.86434196 241.4815 -6431.8631 -2117.6161 -128958699.2767 72647291.8685 31492276.1817 355426.7591 -49491.5168 -31268.686 5.851879255174026 1.7513616764153994 -4.644633372004137 150.0 -1.0
245.148057 40.600859 229.4671 60.5906 109.90413699239146 2460543.98722788 -2490.3566 -6297.04 249.7755 -129118738.9903 72401057.8542 31385551.4949 357166.859 -39415.7836 -25798.8249 5.156047221666995 -2.24657050874157 -5.2211316514345185 150.0 -1.0
245.834503 40.414645 230.5956 60.5906 109.36694032697264 2460543.98952385 -1413.7303 -6580.5999 -783.3318 -129121723.8741 72396454.2992 31383556.1866 357196.0887 -39227.153 -25696.3791 5.655091536009422 -0.6004442841185664 -5.152761917220232 150.0 -1.0
246.063425 40.353658 230.9718 60.5906 109.18754064800652 2460543.99181552 -268.2143 -6533.7727 -1775.1508 -129124702.9662 72391859.2711 31381564.5742 357225.1433 -39038.8636 -25594.117 5.869678894312642 1.072367429278529 -4.825183145781011 150.0 -1.0
246.29239 40.293217 231.348 60.5906 109.00798421442727 2460543.99410718 890.7365 -6159.8221 -2677.8379 -129127681.8665 72387264.1335 31379572.9143 357254.0778 -38850.5616 -25491.8465 5.789979651020758 2.6905526031350746 -4.256189797770276 150.0 -1.0
245.729005 40.630821 230.2941 60.7677 109.56126851046452 2460544.1223988 799.4065 -6175.3227 -2670.8767 -129294139.0024 72129847.1648 31268001.3148 358682.0228 -28291.0172 -19754.248 5.8273127832626415 2.600368197709031 -4.261228152340664 150.0 -1.0
245.482776 40.321351 230.1455 60.4134 109.53088395415786 2460544.24647504 -1425.6521 -6558.6608 -930.1804 -129454552.3146 71880565.7847 31159955.9731 359703.6731 -18051.7448 -14186.0085 5.692577921465727 -0.5100375323854678 -5.1212378192872166 150.0 -1.0
245.711242 40.259838 230.5198 60.4134 109.35312968387876 2460544.2487667 -272.478 -6494.6333 -1912.1339 -129457509.8039 71875958.6232 31157959.1044 359719.2135 -17862.4458 -14083.0225 5.909203734299808 1.1549803757316113 -4.757426227370161 150.0 -1.0
245.670432 40.082999 230.5785 60.2362 109.27580277044919 2460544.25105837 894.3376 -6105.4406 -2798.0917 -129460467.1002 71871351.3539 31155962.1889 359734.633 -17673.1414 -13980.0319 5.82958327867118 2.7612786870304262 -4.155632425449767 150.0 -1.0
244.987301 40.268244 229.4628 60.2362 109.80361294110973 2460544.30169465 -4794.6375 -3830.1624 2870.8175 -129525761.7838 71769522.3935 31111826.7605 360044.5024 -13489.0614 -11703.2976 2.3736968513003376 -6.0423790826234836 -4.091276944259334 150.0 -1.0
245.320872 39.98931 230.1336 60.0591 109.436808306806 2460544.30398632 -4208.7077 -4920.5235 1995.7015 -129528714.6228 71764912.6345 31109828.7664 360057.1303 -13299.6506 -11600.2132 3.521628960779708 -4.926314359472599 -4.713613178448668 150.0 -1.0
245.880258 39.652356 231.1656 59.8819 108.90155622929338 2460544.30627957 -3411.441 -5764.9851 1019.6763 -129531669.3037 71760299.5905 31107829.3484 360069.646 -13110.1053 -11497.0542 4.494268081643895 -3.562747555142778 -5.099225318848539 150.13 -1.0
246.437936 39.320193 232.1872 59.7047 108.37329312748258 2460544.30857397 -2442.555 -6320.5271 -8.1452 -129534625.2822 71755684.1108 31105828.8748 360082.0467 -12920.4607 -11393.8394 5.241975428116348 -2.0188149974092084 -5.22874586693824 150.0 -1.0
246.664572 39.262229 232.5534 59.7047 108.20062261784147 2460544.31086564 -1352.3925 -6558.6542 -1034.3554 -129537577.541 71751074.0295 31103830.741 360094.3117 -12731.0385 -11290.7442 5.725898071800657 -0.3761814325753221 -5.095503403877046 150.0 -1.0
246.788427 39.417264 232.638 59.8819 108.20602365965223 2460544.3131573 -194.4923 -6468.312 -2008.6083 -129540529.6074 71746463.8395 31101832.5601 360106.4558 -12541.6127 -11187.6455 5.923182734482566 1.2854591290033568 -4.706134341454344 150.0 -1.0
247.015567 39.359822 233.0061 59.8819 108.03175717856864 2460544.31544897 973.1513 -6054.1524 -2882.0322 -129543481.4807 71741853.5418 31099834.3325 360118.479 -12352.1833 -11084.5432 5.824024722647891 2.881892888999573 -4.080545976303256 150.0 -1.0
245.442681 40.144213 230.2066 60.2362 109.45190034606792 2460544.56314536 -3466.49 -5751.0355 908.8772 -129861397.5299 71242912.7239 30883579.9802 360705.0358 8129.6749 72.426 4.496169173085036 -3.5216299366038584 -5.1254710591047115 150.2561 -1.0
246.457844 39.688704 231.9836 60.0591 108.56207111863728 2460544.56546631 -2483.7894 -6303.4046 -133.302 -129864365.7715 71238231.6387 30881551.0727 360703.8516 8321.5316 177.0207 5.264934763154088 -1.9639248042494728 -5.22588766470875 150.0 -1.0
247.242746 39.302914 233.3742 59.8819 107.85734021769433 2460544.56777145 -1381.3988 -6531.5956 -1161.7026 -129867313.5787 71233582.3611 30879535.9514 360702.5528 8512.0768 280.9018 5.76063955910513 -0.31757440519868396 -5.059439969766117 150.0 -1.0
246.891246 39.204793 232.9196 59.7047 108.02779950138525 2460544.57006311 -215.9555 -6430.4634 -2125.6841 -129870243.9714 71228960.136 30877532.5554 360701.1404 8701.5052 384.1757 5.964343386639541 1.3354165872806334 -4.638864959210878 150.0 -1.0
245.214971 40.205962 229.8347 60.2362 109.62783754456726 2460544.87579862 -5317.9781 -1641.1353 3862.0071 -130259455.5035 70611343.2048 30609840.8835 359431.2075 33916.1864 14144.7204 0.11987854012659227 -7.12072352989253 -2.855255772028613 150.751 -1.0
245.548187 39.928124 230.5036 60.0591 109.26217747727142 2460544.87809898 -5160.451 -3002.8443 3201.8347 -130262370.8483 70606689.0613 30607823.6533 359413.535 34105.2779 14248.0184 1.4597810570830916 -6.526458629991794 -3.761515849287207 150.1098 -1.0
246.211337 39.378686 231.821 59.7047 108.54580960481452 2460544.88040412 -4741.6405 -4215.4864 2377.9513 -130265292.0547 70602025.1394 30605802.1849 359395.7047 34294.7501 14351.526 2.7294139762398117 -5.60063142483159 -4.479180282592097 150.0 -1.0
246.541269 39.107223 232.4696 59.5276 108.19500787632461 2460544.88269578 -4087.1921 -5209.5139 1439.0632 -130268195.9931 70597388.3695 30603792.4849 359377.8582 34483.1027 14454.4235 3.8552971016315496 -4.399148462689482 -4.966925509736484 150.0 -1.0
245.093598 40.051027 229.7636 60.0591 109.611277981575 2460545.51774351 -5423.4668 -1550.5345 3751.0841 -131065409.6639 69308414.0448 30045116.0435 349857.3162 85986.958 42650.2024 -0.09114383123232983 -7.045669532365082 -3.037580240356874 150.8849 -1.0
245.426407 39.773082 230.4294 59.8819 109.24838253255383 2460545.52004542 -5304.7093 -2900.5405 3057.2276 -131068272.1324 69303727.1201 30043084.6003 349806.4336 86170.35 42750.8161 1.2817643022845004 -6.475515114875634 -3.9117487258545225 150.2654 -1.0
245.758193 39.497268 231.0885 59.7047 108.89042176704355 2460545.52234016 -4919.7348 -4100.6582 2211.4364 -131071125.4893 69299054.688 30041059.4384 349755.5935 86353.1414 42851.1019 2.586786459896757 -5.581464371307447 -4.586344623371588 150.024 -1.0
246.315109 39.165172 232.1053 59.5276 108.36624363063399 2460545.52463318 -4288.3958 -5094.3012 1255.0809 -131073976.5121 69294385.6528 30039035.7489 349704.6759 86535.7666 42951.2981 3.7616890547633375 -4.408687980303512 -5.029638788637385 150.1725 -1.0
246.767529 39.049783 232.834 59.5276 108.02357179290075 2460545.59076572 -3606.2975 -5725.0354 375.9416 -131156118.0445 69159681.8843 29980651.3809 348186.558 91789.958 45834.6522 4.60580966541329 -3.2443133417109493 -5.209763007941888 150.6208 -1.0
246.334264 39.533747 231.9018 59.8819 108.55409946950297 2460566.3306377 2748.0462 -6169.1914 -453.1071 -148333535.8408 23512159.3895 10192646.5869 60450.6968 -330404.1059 -181441.0878 5.304476014944157 1.9944542101256957 5.183448747851571 150.0 -1.0
246.107241 39.592786 231.5337 59.8819 108.72790594380848 2460566.33292936 3720.1722 -5622.7786 575.6639 -148334441.6054 23506847.7049 10190344.0798 60648.6449 -330362.5504 -181420.6678 4.474975391594853 3.503532990216627 5.166345076261223 150.0 -1.0
245.984715 39.437722 231.4547 59.7047 108.71821767113002 2460566.33522103 4505.2041 -4793.5618 1575.4077 -148335347.1441 23501535.9846 10188041.5573 60846.5764 -330320.9004 -181400.196 3.421358375098862 4.838903706768178 4.890442776170403 150.0 -1.0
245.775541 39.867469 230.8736 60.0591 109.0873869178829 2460566.3375127 5063.4761 -3723.1296 2495.6173 -148336252.4568 23496224.2295 10185739.0199 61044.491 -330279.1562 -181379.6723 2.195034685837527 5.931296102385204 4.367201964323613 150.0 -1.0
245.653313 39.712455 230.7975 59.8819 109.0750489004946 2460566.33980436 5366.8422 -2465.2461 3289.8398 -148337157.5434 23490912.4397 10183436.4676 61242.3887 -330237.3175 -181359.0968 0.85626144966227 6.723957288540987 3.622601404067609 150.0 -1.0
246.23037 39.747758 231.6136 60.0591 108.73733231210605 2461241.50566954 80.8637 -5875.336 3310.2949 -68926089.6736 124313385.5569 53887728.8287 -381644.5656 -52712.6952 -47657.4273 6.4674238222374925 2.124309721548343 3.5921801333006442 150.0 -1.0
246.019586 41.116759 230.3846 61.2992 109.67224565440591 2461242.92295722 5011.6638 -1528.656 4244.7671 -72097985.5766 122772149.5903 53219549.2565 -348759.5141 -156965.9018 -101267.8089 3.6289451503398595 6.507881536438163 -1.926901116476074 150.0 -1.0
245.958346 40.569257 230.6723 60.7677 109.38048260792065 2461249.82627521 4949.6982 -3231.4987 3240.9426 -86926835.194 114291297.821 49543447.959 141338.9572 -336080.5202 -170864.7696 5.223685743785609 4.2875306805945 -3.6904627025842274 150.0 -1.0
245.531709 39.557339 230.7223 59.7047 109.06246755771375 2461269.23952818 -2880.2556 -6125.4349 249.9113 -121785716.3128 82645712.9448 35825859.944 -373210.5725 -79749.6564 -62179.5312 5.0079605261908 -2.5665822914015735 -5.222145167672636 150.0 -1.0
246.993827 38.992869 233.1984 59.5276 107.85198378793002 2461269.24186952 -1800.9207 -6479.7992 -803.565 -121789201.509 82641286.3543 35823940.9336 -373170.2549 -79928.532 -62272.45 5.618691781660345 -0.9214528220491047 -5.149608156367549 150.0 -1.0
247.335146 40.211508 232.9412 60.7677 108.29265726866885 2461270.07893123 -769.8714 -6380.5736 -2138.3472 -123022891.3801 81050512.969 35134317.431 -352129.5554 -142009.5116 -94147.8612 6.066087173325059 0.8227227995757858 -4.631040712938557 150.0 -1.0
245.184613 39.462822 230.2834 59.5276 109.22026075444091 2461270.20507846 -1953.8643 -6366.137 -1238.8438 -123206672.937 80809372.4243 35029780.6942 -347852.7714 -150972.8417 -98684.4547 5.741016236315772 -0.779832849284503 -5.035967412190985 150.0 -1.0
245.853003 40.785471 230.3696 60.9449 109.57619911991901 2461418.7407085 1523.7184 -5983.8528 2721.121 56123849.382 -124790935.011 -54094288.6236 390383.5074 -60535.4399 -5777.8035 5.5357508664110355 3.3165895462389017 4.182885463355888 150.5707 -1.0
246.646515 40.387901 231.8067 60.7677 108.83728244843489 2461543.54874711 -2388.6956 -6083.2713 1709.1037 83446631.3446 115788374.8393 50190772.3268 -340152.8286 -151368.229 -100811.1025 5.903110698322962 -0.9523199922568892 4.832041089318864 150.0 -1.0
247.302399 40.967241 232.3953 61.4764 108.74937793880008 2461595.35190587 -114.4507 -6604.0959 -1478.2315 -42127365.5462 134087974.2086 58124125.4211 -365873.3929 59047.5571 1704.9782 5.770283021212208 -1.193038167521081 4.9262879655475444 150.0 -1.0
245.305262 39.617934 230.3561 59.7047 109.23435362008267 2461596.1256083 1664.7978 -6526.7249 628.5529 -44004869.2267 133576770.2274 57902478.0337 -373528.3619 -2948.5039 -30349.0838 5.380243318206794 1.8788754695651733 5.1549264592642885 150.0 -1.0
245.636663 39.342208 231.0121 59.5276 108.87915765680961 2461598.74598264 -3170.1243 -5540.548 -2258.8813 -50305530.828 131678841.4143 57079688.8042 -308469.2095 -200405.5657 -125538.9643 4.4081954003740975 -4.356597078140112 4.525232810415218 150.0 -1.0
245.078792 39.679069 229.9898 59.7047 109.40612540316711 2461600.4335385 4854.9011 -1110.85 4567.3288 -54311775.282 130321907.6391 56491523.2065 -206128.2363 -294768.7233 -165515.1436 1.9699089197547437 7.424480514030165 -0.2871649019743927 150.0 -1.0
246.712425 40.932752 231.5385 61.2992 109.11487982828778 2461602.40634996 -1644.0534 -6383.4906 1519.7525 -58939324.5574 128603758.5468 55746863.6726 -51409.2043 -353904.482 -183540.719 5.900989745675815 -0.3387750005822454 4.908605190053095 150.0 -1.0
245.894591 40.962275 230.3082 61.122 109.65713330366285 2461602.92184524 -378.8109 -6179.425 2720.0475 -60138022.8186 128131612.1876 55542241.8887 -8153.0911 -359150.9786 -182847.9418 6.275354018571306 1.469428603295455 4.184492605568266 150.0 -1.0
246.187671 40.508259 231.0504 60.7677 109.1995872247299 2461602.98570653 -526.2456 -6199.6502 2648.6722 -60286209.699 128072456.0847 55516604.6716 -2774.2932 -359498.5975 -182606.9917 6.274831171352844 1.2925904281473777 4.243323766240792 150.0 -1.0
246.002935 39.807347 231.2436 60.0591 108.91243805412303 2461603.05125107 243.1773 -5969.9369 3163.8344 -60438230.3836 128011588.5416 55490225.7951 2747.3336 -359786.015 -182324.3756 6.309774913703734 2.2613074201944983 3.759058773302968 150.0 -1.0
245.290057 39.246988 230.574 59.3504 109.03689222821579 2461605.56148132 4987.0377 -1685.0858 4236.3263 -66203669.4128 125565168.2511 54430018.032 204273.7056 -319209.5482 -146069.92 3.8310999765326414 6.3658881709307025 -1.9728031843442295 150.0 -1.0
246.088923 39.223662 231.7409 59.5276 108.53737164023607 2461605.68976032 5023.7706 -1653.0778 4205.0654 -66495239.0513 125434153.3093 54373240.5987 213498.5396 -314571.0786 -142997.7228 3.8152046496361436 6.349037222082643 -2.058016316377689 150.0 -1.0
246.644278 38.894615 232.7492 59.3504 108.01907369817943 2461607.42021892 5125.4697 -1922.0498 3961.1873 -70397971.3024 123610410.4391 53582887.7271 320050.318 -231217.2761 -92388.3575 4.225226337521235 5.865681829166078 -2.61568428708415 150.0 -1.0
247.437748 39.999279 233.2288 60.5906 108.1079967796483 2461657.34719441 2695.1079 -5792.4401 2173.8019 -145965420.6817 34437258.8032 14928360.4354 -16319.3402 -356197.5197 -181496.0182 4.756093059577324 3.9473321914465487 4.576954499393651 150.0 -1.0
247.344772 39.091498 233.6521 59.7047 107.68165385763929 2461665.48727145 5046.3547 -1512.7885 4252.2265 -149422081.108 15626487.7885 6775053.3799 393805.429 21055.8484 40815.8131 3.5955344579815574 6.497317850005402 -1.943505208764691 150.0 -1.0
246.082733 40.723897 230.7497 60.9449 109.39390640337118 2461811.5829248 -2778.4889 -6084.1997 930.2737 113811061.6949 -86222560.9045 -37376316.846 -242172.8958 243344.9823 98792.8552 4.866483137743671 -3.023373323840133 -5.129947217242261 150.2429 -1.0
245.410651 39.402247 230.6478 59.5276 109.04976600873304 2461811.77483772 -2786.1193 -6097.6538 810.2227 114132792.18 -85871042.6498 -37223965.6 -255034.4168 232176.8332 92170.1059 4.89196986382467 -2.937139623627445 -5.157306470549656 150.216 -1.0
247.722041 39.556312 233.9256 60.2362 107.68394306247586 2461811.96660447 -2853.2931 -6073.9439 753.2072 114452991.5869 -85518825.6858 -37071311.202 -267271.5656 220456.2281 85329.9642 4.874005302908699 -2.94740901296791 -5.168664981088114 150.2831 -1.0
247.082622 39.902238 232.7651 60.4134 108.28355266899592 2461896.49598071 -963.3688 -6196.5208 -2566.5773 107366918.1802 97257505.9811 42157713.848 -363008.0727 50845.0533 -7617.1076 5.756533610555152 -2.6751263605729036 4.310688128180118 150.0 -1.0
246.417101 40.447796 231.4286 60.7677 109.0184881048976 2461899.58242323 821.1814 -6668.8938 837.412 101771461.8001 102332758.7179 44357687.681 -291673.1268 -194831.5858 -119072.1314 5.570765343018138 1.3352743879477202 5.113377406684094 150.0 -1.0
247.666974 39.942139 233.605 60.5906 107.92753000123446 2461905.60659372 -4377.3321 -5150.4274 -407.4384 90065842.8013 111423475.7327 48298809.5973 184263.2005 -323880.3185 -138400.8692 4.110085855372745 -3.8949299324303865 5.183705589893179 150.0 -1.0
247.564738 40.153831 233.3193 60.7677 108.11084826112875 2461957.83375749 737.9661 -6683.2247 -658.6572 -35326595.6699 135729785.4841 58835510.2332 -22467.4642 -349795.9207 -168840.7383 5.701927409962764 0.1298531303564232 5.154112737973264 150.0 -1.0
245.515711 39.186953 230.9366 59.3504 108.86761911481969 2461957.89911466 1399.9822 -6608.4571 -6.8832 -35487379.114 135694423.3234 58820188.7666 -16782.5877 -350469.2479 -168626.109 5.531801721440388 1.1750367048443877 5.208394729503283 150.0 -1.0
246.31257 40.662861 231.1299 60.9449 109.21140929540945 2461958.78284326 -3960.0001 -3726.9035 -4017.465 -37657148.4979 135200357.1378 58606123.8147 59970.9227 -352168.201 -162202.9683 3.4776609985829054 -6.373973081808799 2.4927126986995063 150.0 -1.0
246.250464 41.054864 230.7692 61.2992 109.48662804384716 2461959.43483687 1347.1283 -6550.3127 941.0329 -39252686.0655 134816887.373 58439976.4179 115305.4416 -344717.0338 -153434.546 5.4497672327397595 1.8620155068128952 5.096602184468111 150.0 -1.0
246.418504 38.95204 232.3867 59.3504 108.1890842038321 2461959.93652816 -3314.4203 -5051.1583 -3032.1521 -40477239.0257 134510888.09 58307393.5108 156324.0444 -334136.046 -144522.6834 4.216848033211136 -5.099575949661812 3.899330169779205 150.0 -1.0
245.74134 39.127451 231.2991 59.3504 108.69823302047337 2461960.00055384 -3300.4563 -5086.409 -2987.8676 -40633313.0207 134471153.7064 58290177.3322 161432.4836 -332492.176 -143257.2799 4.231211532414818 -5.053713887618869 3.9430983798229136 150.0 -1.0
246.48139 40.993533 231.1538 61.2992 109.30085507920921 2461960.95975612 -3523.6076 -5095.2176 -2703.9332 -42965876.3823 133857406.8935 58024247.5324 233636.1292 -300335.9563 -121117.3797 4.033294810187726 -5.0098050846190745 4.200978689682415 150.0 -1.0
245.967004 39.068467 231.6616 59.3504 108.52868870359053 2461961.67202793 258.8961 -6610.7974 1352.8773 -44690869.4311 133379313.9965 57817088.0375 280862.7683 -268027.3523 -101259.1454 5.728600171161133 1.2516035100816225 4.973460048538527 150.0 -1.0
245.862775 39.282674 231.3765 59.5276 108.70834345657407 2461961.7912584 -3647.0213 -5181.3371 -2353.0098 -44979014.3428 133297427.7404 57781605.7093 288146.7023 -261984.5661 -97691.0248 3.9793171426299523 -4.819855744763117 4.464712687431529 150.0 -1.0
246.125044 40.900398 230.6908 61.122 109.47309187166248 2461961.85241015 -4471.1884 -3840.08 -3311.1891 -45126731.1923 133255223.1018 57763317.8338 291808.8972 -258818.8098 -95836.2207 2.6955617094972695 -6.22762067304393 3.592697494409549 150.0 -1.0
246.355484 40.839096 231.0733 61.122 109.28894225224086 2461961.91644618 -4462.0903 -3889.8199 -3264.8553 -45281364.985 133210877.9948 57744102.4012 295589.4648 -255456.2655 -93876.422 2.727971425420037 -6.181705520612013 3.6470232867532104 150.0 -1.0
246.586032 40.77834 231.4559 61.122 109.10458959046838 2461962.04413763 -4532.6483 -3784.651 -3291.1572 -45589559.2301 133121993.9205 57705587.3414 302959.5246 -248608.9116 -89916.5412 2.603103836063383 -6.2531402451168105 3.6164912330968946 150.0 -1.0
246.126057 39.96218 231.3223 60.2362 108.92313219026747 2461962.63277483 825.5874 -6272.9571 2354.5096 -47007593.4731 132704382.8536 57524625.4811 333900.029 -214742.734 -70851.0952 5.767102121208154 2.4436210540716274 4.462294250715786 150.0 -1.0
246.396711 40.078597 231.6424 60.4134 108.81906113254608 2461963.84700222 259.3175 -6254.4347 2527.4227 -49918245.567 131802204.2046 57133667.1355 380528.931 -134966.7945 -28416.3176 6.021747071869754 2.0131588530057343 4.337991907159449 150.0 -1.0
246.943449 40.872554 231.9231 61.2992 108.92880040922387 2461968.12876195 -1172.0599 -5757.6869 3320.5816 -60011189.3534 128187092.6554 55566782.5454 336032.0637 176515.8169 115887.5462 6.769633364535821 0.6982789318408322 3.5808791901808945 150.0 -1.0
246.853963 39.960473 232.3909 60.4134 108.4621919061398 2461968.26853717 5044.9604 -1666.7097 4156.14 -60335789.4756 128057798.9753 55510735.8064 329119.236 185513.9505 119524.896 3.868107778707781 6.286005516581099 -2.1698325132274476 150.0 -1.0
247.232254 40.424316 232.6505 60.9449 108.47998131476675 2461968.77999749 5142.2287 -1639.4297 4046.1101 -61520740.8474 127578668.63 55303036.6762 301192.8988 216874.6788 131840.0113 3.883848272062096 6.181727858333943 -2.427029845607649 150.0 -1.0
246.542392 40.602397 231.51 60.9449 109.02880522590246 2461969.27889607 -1492.8636 -5570.7808 3504.9733 -62672265.6826 127102211.2637 55096490.6547 270181.6233 244779.0986 142215.9977 6.918143732200406 0.26348434721512604 3.347216960219157 150.0 -1.0
246.607899 41.14788 231.2352 61.4764 109.31251081523813 2461969.48250457 5022.6411 -1991.7825 4038.4218 -63140974.777 126905187.0903 55011078.4229 256526.7943 255312.1384 145948.5925 4.305012264512698 5.888079023263799 -2.4457485425227787 150.0 -1.0
246.521339 40.23334 231.7241 60.5906 108.82832024886652 2461969.93014068 5161.1423 -1913.186 3898.9906 -64168859.123 126466792.8951 54821026.5437 224628.8476 276553.312 153064.0931 4.221967064113414 5.825841918933747 -2.7252320431307453 150.0 -1.0
246.839349 41.087095 231.6219 61.4764 109.12495090482 2461970.95235458 5137.4841 -2208.3083 3772.0061 -66502536.0965 125438807.8162 54375362.0164 143490.4488 313978.6068 163247.0488 4.525264345111282 5.477807081865649 -2.951414127633143 150.0 -1.0
245.845688 38.913179 231.5832 59.1732 108.5197426516053 2461971.07986267 5027.2302 -2369.9718 3822.323 -66792280.2661 125307967.0929 54318637.1196 132703.8129 317462.3049 163889.6588 4.711795763932969 5.364353510867659 -2.866017980441029 150.0 -1.0
246.772321 40.542473 231.8902 60.9449 108.8459995436438 2462030.29402331 5197.0175 -1859.6237 3901.7679 -149308342.7923 16681434.4854 7231824.4281 -239574.8431 255964.5898 98753.9761 4.144312310500027 5.861287509551403 -2.7351865145055854 150.0 -1.0
247.047209 40.658529 232.221 61.122 108.7354735308835 2462030.80624205 5244.0496 -1932.1058 3801.9893 -149440899.2725 15484219.5773 6712812.6574 -271639.5091 225034.1095 81226.8271 4.215009001491126 5.724789276980573 -2.9124319807897336 150.0 -1.0
247.997741 39.674763 234.262 60.4134 107.56745835107756 2462031.31842367 5283.7344 -2020.1733 3700.7797 -149561972.8395 14285862.8381 6193309.1777 -299228.7445 190353.9128 62347.4756 4.285410045258686 5.58392762352014 -3.078362375213146 150.0 -1.0
246.876034 40.328544 232.1849 60.7677 108.6558758373069 2462031.57440339 5262.7436 -2116.4254 3676.4321 -149618177.0556 13686514.344 5933485.1844 -311170.9682 171777.2025 52496.7593 4.378710896353368 5.489230202088038 -3.116244377106341 150.0 -1.0
247.311323 39.84455 233.1393 60.4134 108.1047630499479 2462174.23879932 -2869.0819 -6069.0656 793.1017 109331229.9378 -90894868.4058 -39401578.0123 -197069.4369 -283994.9471 -149367.3448 4.8587505300793055 -2.982423282956006 -5.159325211233181 150.2773 -1.0
246.192767 39.009986 232.0242 59.3504 108.35894074841596 2462174.43108379 -2909.3605 -6060.341 707.8512 109671118.7476 -90556277.7514 -39254795.3362 -182616.9995 -293032.3995 -152098.5356 4.860413613728269 -2.9492565987311536 -5.174892646189164 150.2637 -1.0
247.896243 39.885553 233.9812 60.5906 107.74691651250163 2462174.87915376 -3244.1133 -5879.3537 770.1803 110458300.8197 -89763249.2637 -38911006.7054 -147509.7388 -311705.2917 -157230.4518 4.6679840427714225 -3.263689292152776 -5.162518025827012 150.3422 -1.0
247.070849 41.026882 232.0086 61.4764 108.93723946870263 2462175.19888911 -3594.9744 -5643.6638 951.1984 111015855.1645 -89193930.9172 -38664197.7124 -121415.9716 -322919.3201 -159824.4286 4.39793623730111 -3.6765481198744965 -5.123315670992952 150.3534 -1.0
247.768852 39.730819 233.8877 60.4134 107.74673841160009 2462253.15152102 -1211.1436 -5049.711 -4336.9329 121362268.2401 81736272.7992 35429318.4281 -351193.0197 -48435.5681 -55431.467 6.823528124857603 -3.08510699512226 1.6831110300699947 150.0 -1.0
247.118021 39.147872 233.2859 59.7047 107.85477794891013 2462262.69777181 -347.402 -6744.0498 276.7957 105695760.4421 98834877.4684 42841973.2401 322897.4677 -221860.0347 -66968.7278 5.6633989486050496 -0.07635456612661073 5.195740062693517 150.0 -1.0
247.208626 40.056956 232.8527 60.5906 108.28826749754191 2462319.10684862 247.406 -6564.4141 -1642.4043 -25441545.0952 137577555.6066 59636758.8963 384901.2544 -102871.5324 -6666.4917 5.861530949557903 -0.9786556003473379 4.855618813403453 150.0 -1.0
246.979487 40.1152 232.4765 60.5906 108.4684365812141 2462319.74841245 190.256 -6638.4103 -1318.2807 -27040967.7247 137322355.3345 59526164.8031 396522.101 -54823.6005 15723.717 5.783374763091729 -0.8079256480596086 4.980751912410142 150.0 -1.0
247.105538 40.269756 232.563 60.7677 108.47436563632616 2462319.812486 142.5629 -6639.2931 -1321.0232 -27200539.4907 137295994.4829 59514740.4065 397261.822 -49952.9422 17949.2628 5.780481847124008 -0.8511806083515694 4.9824665871574325 150.0 -1.0
246.625286 40.01927 232.0166 60.4134 108.64072700435482 2462319.8788307 1216.2664 -6650.1126 -321.9542 -27365735.2681 137268532.1085 59502838.5136 397946.943 -44899.8399 20249.8788 5.60570354570018 0.7905537790156814 5.188217678746696 150.0 -1.0
247.69237 40.308398 233.4108 60.9449 108.11336607869262 2462319.94056896 16.6658 -6634.6282 -1353.4861 -27519431.6071 137242823.7723 59491696.6991 398510.6085 -40189.534 22386.8685 5.767847519702585 -0.9835783188282221 4.970527282287094 150.0 -1.0
247.002234 40.483123 232.2703 60.9449 108.66308984313291 2462320.0068373 1054.6122 -6674.3801 -389.5385 -27684373.5998 137215065.2539 59479666.2403 399036.3599 -35125.9511 24676.0543 5.633374143335459 0.6039831150883401 5.181799366346809 150.0 -1.0
247.469964 39.246542 233.7423 59.8819 107.68277420634686 2462320.32672669 586.267 -6715.7512 -613.3686 -28480107.5414 137078684.9983 59420557.8573 400420.277 -10604.437 35644.8427 5.692243830799271 0.04264450862137892 5.155318912658231 150.0 -1.0
246.353932 39.902577 231.6942 60.2362 108.74656203862885 2462320.39207401 1159.9081 -6667.0986 -44.6746 -28642563.4754 137050339.2478 59408272.2453 400468.0681 -5585.8416 37866.0031 5.5718285265538725 0.9513103997051076 5.20181847474241 150.0 -1.0
246.816627 40.718145 231.8385 61.122 108.92008344366282 2462320.45908824 2465.0624 -6171.1226 1258.582 -28809128.6764 137021099.2945 59395598.9577 400434.4571 -438.307 40135.9054 4.987255645787627 3.029769582044537 5.000514050686888 150.0 -1.0
246.750391 40.173995 232.1003 60.5906 108.6484546432903 2462320.71351922 1445.7857 -6595.7929 432.9487 -29441197.5337 136908507.2897 59346797.9085 399546.8932 19095.1647 48672.9892 5.457516486480927 1.5529534998040617 5.177141500764151 150.0 -1.0
247.17452 40.812923 232.3077 61.2992 108.74256993242068 2462320.77651841 885.6544 -6709.5277 -50.4945 -29597622.7349 136880242.9655 59334546.9846 399141.8506 23924.8786 50765.1144 5.607114462145546 0.7176541563603678 5.201735673724102 150.0 -1.0
247.277898 40.599462 232.6036 61.122 108.55066476038463 2462320.84177882 1407.5329 -6600.9668 480.4151 -29759628.5206 136850802.7953 59321786.3003 398644.9898 28923.1824 52922.4441 5.463236461835095 1.5585713168184157 5.173147883371598 150.0 -1.0
247.405641 40.75386 232.6923 61.2992 108.55618983272879 2462320.90477846 848.2342 -6714.6948 -2.1897 -29915989.122 136822226.7464 59309400.0684 398090.8867 33742.8437 54995.1144 5.604772019247548 0.7230969162446854 5.201014303481483 150.0 -1.0
247.462259 40.366084 233.0306 60.9449 108.29677171487835 2462321.09496945 -280.8852 -6703.5008 -910.3315 -30387831.3563 136735029.925 59271604.2093 395976.0197 48251.4717 61189.4642 5.667605424368389 -0.9135896564293944 5.097389409246489 150.0 -1.0
247.636809 40.695368 233.0769 61.2992 108.36966154651509 2462321.61011072 604.6131 -6736.6975 228.7071 -31664309.9512 136491855.5275 59166194.9395 386948.344 87075.9225 77426.019 5.6153257755443065 0.6971241051178367 5.194694090523342 150.0 -1.0
247.794435 40.096696 233.6975 60.7677 107.92884388133275 2462321.67232451 -346.8075 -6734.353 -604.226 -31818316.3696 136461795.9207 59153164.5428 385536.0825 91702.4384 79326.9631 5.642762956610454 -0.7369361795463192 5.156186184104854 150.0 -1.0
246.070956 38.854726 231.9439 59.1732 108.35141276177265 2462321.86668987 563.741 -6734.8087 357.5065 -32299237.0127 136366925.8837 59112039.2409 380683.7584 106048.357 85173.7723 5.616744443631205 0.7625317362453643 5.185346321812858 150.0 -1.0
247.508634 40.54096 232.9862 61.122 108.36570681440722 2462322.05884512 399.7833 -6748.0594 330.1177 -32774359.3007 136271706.5252 59070761.7239 375237.778 120051.9093 90809.6151 5.632439604359547 0.6044181975539686 5.186737518226323 150.0 -1.0
247.739417 40.483024 233.3688 61.122 108.18060113241472 2462322.44440913 681.9843 -6679.9477 835.7351 -33726693.5119 136076366.639 58986079.7494 362403.0578 147519.2174 101647.3405 5.601859946459626 1.2289689967838058 5.113525943238528 150.0 -1.0
247.595755 39.401506 233.8335 60.0591 107.68353096519836 2462329.04609969 -1988.4399 -6046.3309 2299.6134 -49783826.9801 131852391.0206 57154628.4732 -116549.6763 338574.4681 137859.8862 6.218106696430823 -0.31583298971781915 4.495656837843853 150.0 -1.0
247.922466 40.251288 233.7909 60.9449 107.92986232155825 2462330.26613771 -983.2823 -5771.9076 3387.3469 -52690301.2234 130891957.1828 56738202.8943 -211927.3478 294886.7724 108995.7536 6.775456380065634 0.9177357850330736 3.4974696432527708 150.0 -1.0
248.024116 40.040135 234.0756 60.7677 107.74674181521529 2462330.32825689 -2125.929 -5806.9612 2748.2083 -52837701.0734 130841576.1569 56716359.2306 -216387.3897 291993.5791 107270.7851 6.444069641288479 -0.3713079646695405 4.159348266816778 150.0 -1.0
246.870149 38.837697 233.1118 59.3504 107.84886369049836 2462330.52405023 -108.9998 -5551.7283 3859.5644 -53301910.1884 130681845.0956 56647105.1776 -230138.9105 282473.5668 101687.6712 6.918374020225571 1.8269670143308758 2.7960470906509727 150.0 -1.0
compile StrayLight.cpp and generate .so file:
Linux:
g++ -fPIC -shared StrayLight.cpp ephcom_wrapper.cpp ephcom.c gnulliver.c -o libstraylight.so
Mac:
g++ -fPIC -shared StrayLight.cpp ephcom_wrapper.cpp ephcom.c gnulliver.c -o libstraylight.dylib
import ctypes
import numpy as np
import astropy.constants as cons
from scipy import interpolate
import math
from astropy.table import Table
import astropy.coordinates as coord
from astropy import units as u
filterPivotWave = {'nuv':2875.5,'u':3629.6,'g':4808.4,'r':6178.2, 'i':7609.0, 'z':9012.9,'y':9627.9}
filterIndex = {'nuv':0,'u':1,'g':2,'r':3, 'i':4, 'z':5,'y':6}
filterCCD = {'nuv':'UV0','u':'UV0','g':'Astro_MB','r':'Astro_MB', 'i':'Basic_NIR', 'z':'Basic_NIR','y':'Basic_NIR'}
bandRange = {'nuv':[2504.0,3230.0],'u':[3190.0,4039.0],'g':[3989.0,5498.0],'r':[5438.0,6956.0], 'i':[6886.0,8469.0], 'z':[8379.0,10855.0],'y':[9217.0, 10900.0], 'GU':[2550, 4000],'GV':[4000, 6200],'GI':[6200,10000]}
Instrument_dir = '/Users/zhangxin/Work/SlitlessSim/CSST_SIM/CSST_C6/straylight/straylight/Instrument/'
SpecOrder = ['-2','-1','0','1','2']
def transRaDec2D(ra, dec):
x1 = np.cos(dec / 57.2957795) * np.cos(ra / 57.2957795);
y1 = np.cos(dec / 57.2957795) * np.sin(ra / 57.2957795);
z1 = np.sin(dec / 57.2957795);
return np.array([x1, y1, z1])
def getAngle132(x1 = 0, y1 = 0, z1 = 0, x2 = 0, y2 = 0, z2 = 0, x3 = 0, y3 = 0, z3 = 0):
cosValue = 0;
angle = 0;
x11 = x1-x3;
y11 = y1-y3;
z11 = z1-z3;
x22 = x2-x3;
y22 = y2-y3;
z22 = z2-z3;
tt = np.sqrt((x11*x11 + y11*y11 + z11* z11) * (x22*x22 + y22*y22 + z22*z22));
if(tt==0):
return 0;
cosValue = (x11*x22+y11*y22+z11*z22)/tt;
if (cosValue > 1):
cosValue = 1;
if (cosValue < -1):
cosValue = -1;
angle = math.acos(cosValue);
return angle * 360 / (2 * math.pi);
def calculateAnglePwithEarth(sat = np.array([0,0,0]), pointing = np.array([0,0,0]), sun = np.array([0,0,0])):
modSat = np.sqrt(sat[0]*sat[0] + sat[1]*sat[1]+sat[2]*sat[2])
modPoint = np.sqrt(pointing[0]*pointing[0] + pointing[1]*pointing[1] + pointing[2]*pointing[2])
withLocalZenithAngle = (pointing[0] * sat[0] + pointing[1] * sat[1] + pointing[2] * sat[2]) / (modPoint*modSat)
innerM_sat_sun = sat[0] * sun[0] + sat[1] * sun[1] + sat[2] * sun[2]
cosAngle = innerM_sat_sun / (modSat * cons.au.value/1000)
isInSunSide = 1
if (cosAngle < -0.3385737): #cos109.79
isInSunSide = -1;
elif cosAngle >= -0.3385737 and cosAngle <= 0.3385737:
isInSunSide = 0;
return math.acos(withLocalZenithAngle)*180/math.pi,isInSunSide
# /**
# * *eCoor = ra, *eCoor+1 = dec
# */
def Cartesian2Equatorial(carCoor = np.array([0,0,0])):
eCoor = np.zeros(2)
if (carCoor[0] > 0 and carCoor[1] >= 0):
eCoor[0] = math.atan(carCoor[1] / carCoor[0]) * 360 / (2 * math.pi)
elif (carCoor[0] < 0):
eCoor[0] = (math.atan(carCoor[1] / carCoor[0]) + math.pi) * 360 / (2 * math.pi)
elif (carCoor[0] > 0 and carCoor[1] < 0):
eCoor[0] = (math.atan(carCoor[1] / carCoor[0]) + 2 * math.pi) * 360 / (2 * math.pi)
elif (carCoor[0] == 0 and carCoor[1] < 0):
eCoor[0] = 270
elif (carCoor[0] == 0 and carCoor[1] > 0):
eCoor[0] = 90
eCoor[1] = math.atan(carCoor[2] / np.sqrt(carCoor[0] * carCoor[0] + carCoor[1] * carCoor[1])) * 360 / (2 * math.pi)
return eCoor
class StrayLight(object):
def __init__(self, jtime = 2460843., sat = np.array([0,0,0]), radec = np.array([0,0])):
self.jtime = jtime
self.sat = sat
self.equator = coord.SkyCoord(radec[0]*u.degree, radec[1]*u.degree,frame='icrs')
self.ecliptic = self.equator.transform_to('barycentrictrueecliptic')
self.pointing = transRaDec2D(radec[0], radec[1])
self.slcdll=ctypes.CDLL('./libstraylight.dylib')
self.slcdll.Calculate.argtypes =[ctypes.c_double ,ctypes.POINTER(ctypes.c_double), ctypes.POINTER(ctypes.c_double),ctypes.POINTER(ctypes.c_double), ctypes.POINTER(ctypes.c_double)]
self.slcdll.PointSource.argtypes =[ctypes.c_double ,ctypes.POINTER(ctypes.c_double), ctypes.POINTER(ctypes.c_double),ctypes.POINTER(ctypes.c_double), ctypes.POINTER(ctypes.c_double)]
self.slcdll.EarthShine.argtypes =[ctypes.c_double ,ctypes.POINTER(ctypes.c_double), ctypes.POINTER(ctypes.c_double),ctypes.POINTER(ctypes.c_double), ctypes.POINTER(ctypes.c_double)]
self.slcdll.Zodiacal.argtypes =[ctypes.c_double ,ctypes.POINTER(ctypes.c_double), ctypes.POINTER(ctypes.c_double)]
self.slcdll.ComposeY.argtypes=[ctypes.POINTER(ctypes.c_double),ctypes.POINTER(ctypes.c_double),ctypes.POINTER(ctypes.c_double)]
self.slcdll.Init()
def getFilterAndCCD_Q(self, filter = 'i'):
ccd_fn = Instrument_dir + 'ccd/' + filterCCD[filter] + '.txt'
filter_fn = Instrument_dir + 'filters/' + filter + '.txt'
q_ccd_f = np.loadtxt(ccd_fn)
q_fil_f = np.loadtxt(filter_fn)
band_s = 2000
band_e = 11000
q_ccd = np.zeros([q_ccd_f.shape[0]+2,q_ccd_f.shape[1]])
q_ccd[1:-1,:] = q_ccd_f
q_ccd[0] = [band_s,0]
q_ccd[-1] = [band_e,0]
q_fil = np.zeros([q_fil_f.shape[0]+2,q_fil_f.shape[1]])
q_fil[1:-1,:] = q_fil_f
q_fil[0] = [band_s,0]
q_fil[-1] = [band_e,0]
q_fil_i = interpolate.interp1d(q_fil[:,0], q_fil[:,1])
q_ccd_i = interpolate.interp1d(q_ccd[:,0], q_ccd[:,1])
bands = np.arange(bandRange[filter][0], bandRange[filter][1],0.5)
q_ccd_fil = q_fil_i(bands)*q_ccd_i(bands)
return np.trapz(q_ccd_fil, bands)/(bandRange[filter][1]-bandRange[filter][0])
def caculateEarthShineFilter(self, filter = 'i', pixel_size_phy = 10 ):
sat = (ctypes.c_double*3)()
sat[:] = self.sat
ob = (ctypes.c_double*3)()
ob[:]=self.pointing
py1 = (ctypes.c_double*3)()
py2 = (ctypes.c_double*3)()
self.slcdll.ComposeY(ob,py1,py2)
earth_e1 = (ctypes.c_double*7)()
self.slcdll.EarthShine(self.jtime,sat,ob,py1,earth_e1)
earth_e2 = (ctypes.c_double*7)()
self.slcdll.EarthShine(self.jtime,sat,ob,py2,earth_e2)
band_earth_e1 = earth_e1[:][filterIndex[filter]]
band_earth_e2 = earth_e2[:][filterIndex[filter]]
q=self.getFilterAndCCD_Q(filter=filter)
p_lambda = filterPivotWave[filter]
c = cons.c.value
h = cons.h.value
pix_earth_e1 = band_earth_e1/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
pix_earth_e2 = band_earth_e2/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
if pix_earth_e1< pix_earth_e2:
return pix_earth_e1, py1[:]
else:
return pix_earth_e2, py2[:]
"""
calculate zodiacal call c++ program, seems to have some problem
"""
def calculateZodiacalFilter1(self, filter = 'i', pixel_size_phy = 10 ):
sat = (ctypes.c_double*3)()
sat[:] = self.sat
ob = (ctypes.c_double*3)()
ob[:]=self.pointing
zodical_e = (ctypes.c_double*7)()
self.slcdll.Zodiacal(self.jtime,ob,zodical_e)
ob1 = (ctypes.c_double*2)()
ob1[:] = np.array([self.ecliptic.lon.value, self.ecliptic.lat.value])
zodical_e1 = (ctypes.c_double*7)()
self.slcdll.Zodiacal1(ob1,zodical_e1)
band_zodical_e = zodical_e[:][filterIndex[filter]]
q=self.getFilterAndCCD_Q(filter=filter)
p_lambda = filterPivotWave[filter]
c = cons.c.value
h = cons.h.value
pix_zodical_e = band_zodical_e/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
return pix_zodical_e, band_zodical_e
"""
calculate zodiacal use python
"""
def calculateZaodiacalFilter2(self,filter = 'i', aper = 2, pixelsize = 0.074, sun_pos = np.array([0,0,0])):
spec, v_mag = self.calculateZaodicalSpec(longitude = self.ecliptic.lon.value, latitude = self.ecliptic.lat.value, sun_pos = sun_pos)
# spec = self.calculateZaodicalSpec(longitude = lon, latitude = lat)
throughputFn = Instrument_dir + 'throughputs/' + filter + '_throughput.txt'
throughput = np.loadtxt(throughputFn)
deltL = 0.5
lamb = np.arange(bandRange[filter][0], bandRange[filter][1], deltL)
speci = interpolate.interp1d(spec['WAVELENGTH'], spec['FLUX'])
y = speci(lamb)
# erg/s/cm2/A --> photo/s/m2/A
flux = y * lamb / (cons.h.value * cons.c.value) * 1e-13
throughput_i = interpolate.interp1d(throughput[:, 0], throughput[:, 1])
throughput_ = throughput_i(lamb)
sky_pix = np.trapz(flux*throughput_, lamb) * math.pi * aper*aper/4 * pixelsize * pixelsize
# sky_pix_e = np.trapz(y, lamb) * math.pi * aper*aper/4 * pixelsize * pixelsize/(10*10*1e-6*1e-6)*1e-7*1e4
return sky_pix, v_mag#, sky_pix_e
def caculateStarLightFilter(self, filter = 'i', pointYaxis = np.array([1,1,1]), pixel_size_phy = 10 ):
sat = (ctypes.c_double*3)()
sat[:] = self.sat
ob = (ctypes.c_double*3)()
ob[:]=self.pointing
py = (ctypes.c_double*3)()
py[:] = pointYaxis
q=self.getFilterAndCCD_Q(filter=filter)
p_lambda = filterPivotWave[filter]
c = cons.c.value
h = cons.h.value
star_e1 = (ctypes.c_double*7)()
self.slcdll.PointSource(self.jtime,sat,ob,py,star_e1)
band_star_e1 = star_e1[:][filterIndex[filter]]
pix_star_e1 = band_star_e1/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
return pix_star_e1
def caculateStrayLightFilter(self, filter = 'i', pixel_size_phy = 10 ):
sat = (ctypes.c_double*3)()
sat[:] = self.sat
ob = (ctypes.c_double*3)()
ob[:]=self.pointing
py1 = (ctypes.c_double*3)()
py2 = (ctypes.c_double*3)()
self.slcdll.ComposeY(ob,py1,py2)
earth_e1 = (ctypes.c_double*7)()
self.slcdll.EarthShine(self.jtime,sat,ob,py1,earth_e1)
earth_e2 = (ctypes.c_double*7)()
self.slcdll.EarthShine(self.jtime,sat,ob,py2,earth_e2)
zodical_e = (ctypes.c_double*7)()
self.slcdll.Zodiacal(self.jtime,ob,zodical_e)
band_earth_e1 = earth_e1[:][filterIndex[filter]]
band_earth_e2 = earth_e2[:][filterIndex[filter]]
# band_earth_e1 = 0
# band_earth_e2 = 0
band_zodical_e = zodical_e[:][filterIndex[filter]]
q=self.getFilterAndCCD_Q(filter=filter)
p_lambda = filterPivotWave[filter]
c = cons.c.value
h = cons.h.value
pix_earth_e1 = band_earth_e1/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
pix_earth_e2 = band_earth_e2/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
pix_zodical_e = band_zodical_e/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
# star_e1 = (ctypes.c_double*7)()
# self.slcdll.PointSource(self.jtime,sat,ob,py1,star_e1)
# # star_e2 = (ctypes.c_double*7)()
# # self.slcdll.PointSource(self.jtime,sat,ob,py2,star_e2)
# band_star_e1 = star_e1[:][filterIndex[filter]]
# # band_star_e2 = star_e2[:][filterIndex[filter]]
# pix_star_e1 = band_star_e1/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
pix_star_e1 = 0
# pix_star_e2 = band_star_e2/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
return pix_earth_e1+pix_zodical_e+pix_star_e1, pix_zodical_e, pix_earth_e1, pix_earth_e2
def caculateStrayLightGrating(self, grating = 'GU', pixel_size_phy = 10, normFilter = 'g', aper = 2, pixelsize = 0.074):
sat = (ctypes.c_double*3)()
sat[:] = self.sat
ob = (ctypes.c_double*3)()
ob[:]=self.pointing
py1 = (ctypes.c_double*3)()
py2 = (ctypes.c_double*3)()
self.slcdll.ComposeY(ob,py1,py2)
earth_e1 = (ctypes.c_double*7)()
self.slcdll.EarthShine(self.jtime,sat,ob,py1,earth_e1)
earth_e2 = (ctypes.c_double*7)()
self.slcdll.EarthShine(self.jtime,sat,ob,py2,earth_e2)
# zodical_e = (ctypes.c_double*7)()
# self.slcdll.Zodiacal(self.jtime,ob,zodical_e)
band_earth_e1 = earth_e1[:][filterIndex[normFilter]]
band_earth_e2 = earth_e2[:][filterIndex[normFilter]]
band_earth_e = np.min([band_earth_e1, band_earth_e2])
# band_earth_e1 = 0
# band_earth_e2 = 0
# band_zodical_e = zodical_e[:][filterIndex[normFilter]]
q=self.getFilterAndCCD_Q(filter=normFilter)
p_lambda = filterPivotWave[normFilter]
c = cons.c.value
h = cons.h.value
pix_earth_e = band_earth_e/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
# pix_earth_e2 = band_earth_e2/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
# pix_zodical_e = band_zodical_e/(h*c/(p_lambda*1e-10))*pixel_size_phy*1e-6*pixel_size_phy*1e-6*q
# pix_earth_e = np.min([pix_earth_e1, pix_earth_e2])
# zodical_v, zodical_spec = self.calculatSkylightBySpec(specType = 'zodical', filter = 'g', aper = 2, pixelsize = 0.074)
earthshine_v, earthshine_spec = self.calculatSkylightBySpec(specType = 'earthshine', filter = 'g', aper = aper, pixelsize = pixelsize)
lamb_earth = earthshine_spec['WAVELENGTH']
flux_earth = earthshine_spec['FLUX']*pix_earth_e/earthshine_v
print(pix_earth_e,earthshine_v)
earth_v_grating = 0
for s_order in SpecOrder:
thpFn = Instrument_dir + 'sls_conf/' + grating + '.Throughput.' + s_order + 'st.fits'
thp_ = Table.read(thpFn)
thpFn_i = interpolate.interp1d(thp_['WAVELENGTH'], thp_['SENSITIVITY'])
thp = thpFn_i(lamb_earth)
beamsEarth = np.trapz(flux_earth*thp,lamb_earth)* math.pi*aper*aper/4 * pixelsize * pixelsize
earth_v_grating = earth_v_grating + beamsEarth
print(beamsEarth)
# print(earthshine_v, pix_earth_e, earth_v_grating)
return earth_v_grating
def calculatSkylightBySpec(self, specType = 'earthshine', filter = 'g', aper = 2, pixelsize = 0.074, s = 2000, e = 11000):
specFn = ''
if specType == 'zodical':
specFn=Instrument_dir + 'sky/zodiacal.dat'
elif specType == 'earthshine':
specFn= Instrument_dir + 'sky/earthShine.dat'
spec = np.loadtxt(specFn)
throughputFn = Instrument_dir + 'throughputs/' + filter + '_throughput.txt'
throughput = np.loadtxt(throughputFn)
deltL = 0.5
lamb = np.arange(bandRange[filter][0], bandRange[filter][1], deltL)
speci = interpolate.interp1d(spec[:, 0], spec[:, 1])
y = speci(lamb)
# erg/s/cm2/A --> photo/s/m2/A
flux = y * lamb / (cons.h.value * cons.c.value) * 1e-13
throughput_i = interpolate.interp1d(throughput[:, 0], throughput[:, 1])
throughput_ = throughput_i(lamb)
sky_pix = np.trapz(flux*throughput_, lamb) * math.pi * aper*aper/4 * pixelsize * pixelsize
lamb = np.arange(s, e, deltL)
speci = interpolate.interp1d(spec[:, 0], spec[:, 1])
y = speci(lamb)
# erg/s/cm2/A --> photo/s/m2/A
flux = y * lamb / (cons.h.value * cons.c.value) * 1e-13
return sky_pix, Table(np.array([lamb, flux]).T,names=('WAVELENGTH', 'FLUX'))
def calculateZaodicalSpec(self,longitude = 50, latitude = 60, sun_pos = np.array([0,0,0])):
from scipy.interpolate import interp2d
from scipy.interpolate import griddata
z_map_fn = Instrument_dir + 'Zodiacal_map1.dat'
ZL = np.loadtxt(z_map_fn)
# zl_sh = ZL.shape
# x = np.arange(0,zl_sh[1],1)
# y = np.arange(0,zl_sh[0],1)
x = ZL[0,1:]
y = ZL[1:,0]
X,Y = np.meshgrid(x,y)
# f_sur = interp2d(X,Y,ZL,kind='linear')
sun_radec = Cartesian2Equatorial(sun_pos)
sun_eclip = coord.SkyCoord(sun_radec[0]*u.degree, sun_radec[1]*u.degree,frame='icrs')
sun_equtor = sun_eclip.transform_to('barycentrictrueecliptic')
longitude = longitude - (sun_equtor.lon*u.degree).value
longitude = np.abs(longitude)
print((sun_equtor.lon*u.degree).value)
if (longitude > 180):
longitude = 360 - longitude
latitude = np.abs(latitude)
lo = longitude
la = latitude
zl = griddata((X.flatten(),Y.flatten()),ZL[1:,1:].flatten(),(la,lo), method='cubic').min()
zl = zl*(math.pi*math.pi)/(180*180)/(3600*3600)*1e-4*1e7*1e-8*1e-4
# print(zl , '\n')
zodical_fn = Instrument_dir + 'sky/zodiacal.dat'
spec = np.loadtxt(zodical_fn)
speci = interpolate.interp1d(spec[:, 0], spec[:, 1])
flux5000 = speci(5000)
f_ration = zl/flux5000
v_mag = np.log10(f_ration)*(-2.5)+22.1
# print("factor:", v_mag, lo, la)
return Table(np.array([spec[:,0], spec[:,1]*f_ration]).T,names=('WAVELENGTH', 'FLUX')), v_mag
def testZodiacal(lon = 285.04312526255366, lat = 30.):
c_eclip = coord.SkyCoord(lon*u.degree, lat*u.degree,frame='barycentrictrueecliptic')
c_equtor = c_eclip.transform_to('icrs')
sl = StrayLight(jtime = 2459767.00354975, sat = np.array([]), radec = np.array([(c_equtor.ra*u.degree).value, (c_equtor.dec*u.degree).value]))
e_zol, v_mag = sl.calculateZaodiacalFilter2(filter = 'i', sun_pos=np.array([-3.70939436e+07, 1.35334903e+08, 5.86673104e+07]))
print(e_zol)
# ju=2.4608437604166665e+06
# sat = (ctypes.c_double*3)()
# sat[:] = np.array([5873.752, -1642.066, 2896.744])
# ob = (ctypes.c_double*3)()
# ob[:]=np.array([0.445256,0.76061,-0.47246])
# sl = StrayLight(jtime = ju, sat = np.array([5873.752, -1642.066, 2896.744]), pointing = np.array([-0.445256,-0.76061,0.47246]))
fn = '/Users/zhangxin/Work/SurveyPlan/point/csst_survey_sim_20211028/E17.5_b17.5_beta_11.6_opt_transtime_1_CMG_1_dp_2_0.25_da_10_Texp_1.5_DEC60_500_0.1_800_1000_+5deg.dat'
surveylist = np.loadtxt(fn)
sky_pix = np.zeros([surveylist.shape[0],7])
i = 693438
c_eclip = coord.SkyCoord(surveylist[:,2]*u.degree, surveylist[:,1]*u.degree,frame='barycentrictrueecliptic')
c_equtor = c_eclip.transform_to('icrs')
# pointing = transRaDec2D((c_equtor[i].ra*u.degree).value, (c_equtor[i].dec*u.degree).value)
# # print(ju, pointing, surveylist[i,3:9])
# ju = surveylist[i,0]
# sl = StrayLight(jtime = ju, sat = surveylist[i,3:6], pointing = pointing)
# sl.caculateStrayLightGrating(grating = 'GI', pixel_size_phy = 10, normFilter = 'g')
for i in np.arange(surveylist.shape[0]):
print(i)
if i > 300:
break
# if i != 300:
# continue
# if i != 693438:
# continue
ju = surveylist[i,0]
pointing = transRaDec2D((c_equtor[i].ra*u.degree).value, (c_equtor[i].dec*u.degree).value)
# print(ju, pointing, surveylist[i,3:9])
sl = StrayLight(jtime = ju, sat = surveylist[i,3:6], radec = np.array([(c_equtor[i].ra*u.degree).value, (c_equtor[i].dec*u.degree).value]))
# strayl_i,s_zoldical ,s_earth, s_earth1 = sl.caculateStrayLightFilter(filter = 'i')
# print(i,strayl_i,s_zoldical,s_earth, s_earth1)
p_cart= transRaDec2D((c_equtor[i].ra*u.degree).value, (c_equtor[i].dec*u.degree).value)
sky_pix[i,6] = getAngle132(x1 = surveylist[i,6], y1 = surveylist[i,7], z1 = surveylist[i,8], x2 = p_cart[0], y2 = p_cart[1], z2 = p_cart[2], x3 = 0, y3 = 0, z3 = 0)
earthZenithAngle,isInSunSide = calculateAnglePwithEarth(sat = surveylist[i,3:6], pointing = pointing, sun = surveylist[i,6:9])
sky_pix[i,4] = earthZenithAngle
sky_pix[i,5] = isInSunSide
e1,py = sl.caculateEarthShineFilter(filter = 'i')
# e2, e2_ = sl.calculateZodiacalFilter1(filter = 'i')
e3, v_mag = sl.calculateZaodiacalFilter2(filter = 'i', sun_pos=surveylist[i,6:9])
e4 = sl.caculateStarLightFilter(filter = 'i',pointYaxis = py)
# e4 = 0
# e4 = sl.caculateStrayLightGrating(grating = 'GI', pixel_size_phy = 10, normFilter = 'g')
sky_pix[i,0] = e1
sky_pix[i,1] = e3
sky_pix[i,2] = e4
sky_pix[i,3] = v_mag
print(e1,e3,e4)
# print(e1,e2,e3,e4)
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment