Projects
mass
atmos
Log In
Username
Password
We truncated the diff of some files because they were too big. If you want to see the full diff for every file,
click here
.
Overview
Repositories
Revisions
Requests
Users
Attributes
Meta
Expand all
Collapse all
Changes of Revision 13
View file
atmos.spec
Changed
@@ -1,5 +1,5 @@ # -# spec file for package atmos (Version 2.97.2) +# spec file for package atmos (Version 2.97.3) # # Copyright (c) 2011 Matwey V. Kornilov. # This file and all modifications and additions to the pristine @@ -13,7 +13,7 @@ BuildRequires: make gcc gcc-c++ gsl-devel boost-devel Summary: MASS/DIMM data processing tool Name: atmos -Version: 2.97.2 +Version: 2.97.3 Release: 1 Source: %{name}-%{version}.tar.gz Source1: spectra.tar.gz @@ -60,6 +60,8 @@ %attr(644,root,root) %{_datadir}/atmos/spectra/* %changelog +* Fri Mar 04 2011 - matwey.kornilov@gmail.com +- version 2.97.3 * Fri Feb 04 2011 - matwey.kornilov@gmail.com - version 2.97.2 * Mon Jan 31 2011 - matwey.kornilov@gmail.com
View file
atmos-2.97.2.tar.gz/configure -> atmos-2.97.3.tar.gz/configure
Changed
@@ -1,6 +1,6 @@ #! /bin/sh # Guess values for system-dependent variables and create Makefiles. -# Generated by GNU Autoconf 2.63 for atmos 2.97.2. +# Generated by GNU Autoconf 2.63 for atmos 2.97.3. # # Report bugs to <mass-general@curl.sai.msu.ru>. # @@ -596,8 +596,8 @@ # Identity of this package. PACKAGE_NAME='atmos' PACKAGE_TARNAME='atmos' -PACKAGE_VERSION='2.97.2' -PACKAGE_STRING='atmos 2.97.2' +PACKAGE_VERSION='2.97.3' +PACKAGE_STRING='atmos 2.97.3' PACKAGE_BUGREPORT='mass-general@curl.sai.msu.ru' ac_unique_file="src/atmos.h" @@ -1293,7 +1293,7 @@ # Omit some internal or obsolete options to make the list less imposing. # This message is too long to be a string in the A/UX 3.1 sh. cat <<_ACEOF -\`configure' configures atmos 2.97.2 to adapt to many kinds of systems. +\`configure' configures atmos 2.97.3 to adapt to many kinds of systems. Usage: $0 OPTION... VAR=VALUE... @@ -1359,7 +1359,7 @@ if test -n "$ac_init_help"; then case $ac_init_help in - short | recursive ) echo "Configuration of atmos 2.97.2:";; + short | recursive ) echo "Configuration of atmos 2.97.3:";; esac cat <<\_ACEOF @@ -1449,7 +1449,7 @@ test -n "$ac_init_help" && exit $ac_status if $ac_init_version; then cat <<\_ACEOF -atmos configure 2.97.2 +atmos configure 2.97.3 generated by GNU Autoconf 2.63 Copyright (C) 1992, 1993, 1994, 1995, 1996, 1998, 1999, 2000, 2001, @@ -1463,7 +1463,7 @@ This file contains any messages produced by compilers while running configure, to aid debugging if configure makes a mistake. -It was created by atmos $as_me 2.97.2, which was +It was created by atmos $as_me 2.97.3, which was generated by GNU Autoconf 2.63. Invocation command line was $ $0 $@ @@ -2312,7 +2312,7 @@ # Define the identity of the package. PACKAGE='atmos' - VERSION='2.97.2' + VERSION='2.97.3' cat >>confdefs.h <<_ACEOF @@ -5455,7 +5455,7 @@ # report actual input values of CONFIG_FILES etc. instead of their # values after options handling. ac_log=" -This file was extended by atmos $as_me 2.97.2, which was +This file was extended by atmos $as_me 2.97.3, which was generated by GNU Autoconf 2.63. Invocation command line was CONFIG_FILES = $CONFIG_FILES @@ -5518,7 +5518,7 @@ _ACEOF cat >>$CONFIG_STATUS <<_ACEOF || ac_write_fail=1 ac_cs_version="\\ -atmos config.status 2.97.2 +atmos config.status 2.97.3 configured by $0, generated by GNU Autoconf 2.63, with options \\"`$as_echo "$ac_configure_args" | sed 's/^ //; s/\\""\`\$/\\\\&/g'`\\"
View file
atmos-2.97.2.tar.gz/configure.in -> atmos-2.97.3.tar.gz/configure.in
Changed
@@ -2,7 +2,7 @@ # Process this file with autoconf to produce a configure script. AC_PREREQ(2.59) -AC_INIT(atmos, 2.97.2, mass-general@curl.sai.msu.ru) +AC_INIT(atmos, 2.97.3, mass-general@curl.sai.msu.ru) AM_INIT_AUTOMAKE() AC_CONFIG_HEADERS(src/config.h) AC_CONFIG_SRCDIR(src/atmos.h)
View file
atmos-2.97.2.tar.gz/src/covariance_decomposition.h -> atmos-2.97.3.tar.gz/src/covariance_decomposition.h
Changed
@@ -1,5 +1,5 @@ /* - $Id: covariance_decomposition.h 182 2011-01-31 11:26:03Z matwey $ + $Id: covariance_decomposition.h 206 2011-02-15 09:40:12Z matwey $ Copyright (C) 2011 Sternberg Astronomical Institute, MSU This program is free software: you can redistribute it and/or modify @@ -91,6 +91,7 @@ value_type min = std::abs( weight(0,0) ); value_type max = std::abs( weight(0,0) ); for( typename triangular_matrix< value_type >::size_type i = 1; i < weight.size1(); ++i ) { + if( !std::isfinite(weight(i,i)) ) throw std::logic_error( "Infinite diagonal value" ); if( std::abs( weight(i,i) ) > max ) max = std::abs( weight(i,i) ); if( std::abs( weight(i,i) ) < min ) min = std::abs( weight(i,i) ); }
View file
atmos-2.97.2.tar.gz/src/file_output.cpp -> atmos-2.97.3.tar.gz/src/file_output.cpp
Changed
@@ -1,5 +1,5 @@ /* - $Id: file_output.cpp 200 2011-02-03 14:45:13Z matwey $ + $Id: file_output.cpp 207 2011-03-01 19:02:17Z matwey $ Copyright (C) 2009 Sternberg Astronomical Institute, MSU This program is free software: you can redistribute it and/or modify @@ -45,17 +45,17 @@ boost::io::ios_all_saver asaver( m_stm ); m_stm << "M " << time << " Background measurements" << std::endl; } -void file::atmosphere( const time_type& time, double free_seeing, double free_seeing_err, double dimm_seeing, double dimm_seeing_err, double dimm_turbulence, double dimm_turbulence_err, double isoplanatic_angle, double isoplanatic_angle_err, double altitude_eff, double altitude_eff_err, double desi_time_constant, double desi_time_constant_err, double time_constant, double time_constant_err, double moment_0, double moment_0_err, double moment_2, double moment_2_err ) { +void file::atmosphere( const time_type& time, double free_seeing, double free_seeing_rerr, double dimm_seeing, double dimm_seeing_rerr, double dimm_turbulence, double dimm_turbulence_rerr, double isoplanatic_angle, double isoplanatic_angle_rerr, double altitude_eff, double altitude_eff_rerr, double desi_time_constant, double desi_time_constant_rerr, double time_constant, double time_constant_rerr, double moment_0, double moment_0_rerr, double moment_2, double moment_2_rerr ) { static datum skipped; boost::io::ios_all_saver asaver( m_stm ); m_stm << "A " << time << " " - << std::fixed << datum( free_seeing, free_seeing_err/free_seeing ) << " " << datum( dimm_seeing, dimm_seeing_err/dimm_seeing ) << " " - << std::scientific << datum( moment_0, moment_0_err/moment_0 ) << " " << datum( dimm_turbulence, dimm_turbulence_err/dimm_turbulence ) << " " - << std::fixed << datum( altitude_eff, altitude_eff_err/altitude_eff ) << " " << skipped << " " - << datum( isoplanatic_angle, isoplanatic_angle_err/isoplanatic_angle ) << " " - << std::scientific << datum( moment_2, moment_2_err/moment_2 ) << " " - << std::fixed << datum( desi_time_constant, desi_time_constant_err/time_constant ) << " " - << std::fixed << datum( time_constant, time_constant_err/time_constant ) << std::endl; + << std::fixed << datum( free_seeing, free_seeing_rerr ) << " " << datum( dimm_seeing, dimm_seeing_rerr ) << " " + << std::scientific << datum( moment_0, moment_0_rerr ) << " " << datum( dimm_turbulence, dimm_turbulence_rerr ) << " " + << std::fixed << datum( altitude_eff, altitude_eff_rerr ) << " " << skipped << " " + << datum( isoplanatic_angle, isoplanatic_angle_rerr ) << " " + << std::scientific << datum( moment_2, moment_2_rerr ) << " " + << std::fixed << datum( desi_time_constant, desi_time_constant_rerr ) << " " + << std::fixed << datum( time_constant, time_constant_rerr ) << std::endl; } void file::profile( const time_type& time, const std::string& method, double cn2seeing, double chi2, const vector_type& cn2, const vector_type& ecn2, const vector_type& grid ){ assert( cn2.size() == ecn2.size() );
View file
atmos-2.97.2.tar.gz/src/file_output.h -> atmos-2.97.3.tar.gz/src/file_output.h
Changed
@@ -1,5 +1,5 @@ /* - $Id: file_output.h 200 2011-02-03 14:45:13Z matwey $ + $Id: file_output.h 207 2011-03-01 19:02:17Z matwey $ Copyright (C) 2009 Sternberg Astronomical Institute, MSU This program is free software: you can redistribute it and/or modify @@ -50,7 +50,7 @@ void comment( const std::string& comment ); void normal( const time_type& time, double mz ); void background( const time_type& time ); - void atmosphere( const time_type& time, double free_seeing, double free_seeing_err, double dimm_seeing, double dimm_seeing_err, double dimm_turbulence, double dimm_turbulence_err, double isoplanatic_angle, double isoplanatic_angle_err, double altitude_eff, double altitude_eff_err, double desi_time_constant, double desi_time_constant_err, double time_constant, double time_constant_err, double moment_0, double moment_0_err, double moment_2, double moment_2_err ); + void atmosphere( const time_type& time, double free_seeing, double free_seeing_rerr, double dimm_seeing, double dimm_seeing_rerr, double dimm_turbulence, double dimm_turbulence_rerr, double isoplanatic_angle, double isoplanatic_angle_rerr, double altitude_eff, double altitude_eff_rerr, double desi_time_constant, double desi_time_constant_rerr, double time_constant, double time_constant_rerr, double moment_0, double moment_0_rerr, double moment_2, double moment_2_rerr ); void profile( const time_type& time, const std::string& method, double cn2seeing, double chi2, const vector_type& cn2, const vector_type& ecn2, const vector_type& grid ); void residual( const time_type& time, const std::string& method, const vector_type& o_csi ); void fluxes( const time_type& time, const vector_type& result_flux );
View file
atmos-2.97.2.tar.gz/src/indices.h -> atmos-2.97.3.tar.gz/src/indices.h
Changed
@@ -1,5 +1,5 @@ /* - $Id: indices.h 200 2011-02-03 14:45:13Z matwey $ + $Id: indices.h 205 2011-02-04 11:02:23Z matwey $ Copyright (C) 2009 Sternberg Astronomical Institute, MSU This program is free software: you can redistribute it and/or modify @@ -61,8 +61,8 @@ timestamp_type timestamp; flux_type flux; indices_type indices; - indices_type wind_indices; indesis_type indesis; + indices_type wind_indices; sample( const timestamp_type& time, const flux_type& flux,
View file
atmos-2.97.2.tar.gz/src/integrals.cpp -> atmos-2.97.3.tar.gz/src/integrals.cpp
Changed
@@ -1,5 +1,5 @@ /* - $Id: integrals.cpp 200 2011-02-03 14:45:13Z matwey $ + $Id: integrals.cpp 210 2011-03-02 18:59:06Z matwey $ Copyright (C) 2009 Sternberg Astronomical Institute, MSU This program is free software: you can redistribute it and/or modify @@ -36,8 +36,73 @@ namespace atmos { namespace integrals { +value_type seeing( const value_type& moment0 ) { + const value_type power = 3.0/5.0; + const value_type k_lambda = 1.1486983549970351 * 15.848931924611138; // lambda^-1/5 + const value_type k = 5.307 * 206264.8; + + return k*k_lambda*std::pow( moment0, power ); +} +value_type effective_altitude( const value_type& moment0, const value_type& moment1 ) { + const value_type k = 0.001; // km in meter + + return k*moment1/moment0; +} +value_type isoplanatic_angle( const value_type& moment53 ) { + const value_type power = -3.0/5.0; + const value_type k_lambda = 0.43527528164806206 * 6.3095734448019363e-08; // lambda^6/5 + const value_type k = 0.058 * 206264.8; + + return k*k_lambda*std::pow( moment53, power ); +} +value_type time_constant_desi( const value_type& desi ) { + const value_type power = -3.0/5.0; + const value_type k = 0.175; + + return k*std::pow( desi, power ); +} +value_type time_constant( const value_type& wind_moment0, const value_type& moment0 ) { + const value_type power = -0.1; + const value_type k_lambda = 0.43527528164806206 * 6.3095734448019363e-08; // lambda^6/5 + const value_type k = 0.058 * 1e3; + const value_type k_wind = 6e6 / 6; + + return k*k_lambda*std::pow( moment0, power )/std::sqrt( k_wind*wind_moment0 ); +} + +std::pair< value_type, value_type > seeing( const value_type& moment0, const value_type& moment0_rerr ) { + return std::make_pair( seeing(moment0), 3.0/5.0*moment0_rerr ); +} +std::pair< value_type, value_type > effective_altitude( const value_type& moment0, const value_type& moment0_rerr, const value_type& moment1, const value_type& moment1_rerr ) { + return std::make_pair( effective_altitude(moment0,moment1), std::sqrt(moment0_rerr*moment0_rerr+moment1_rerr*moment1_rerr) ); +} +std::pair< value_type, value_type > isoplanatic_angle( const value_type& moment53, const value_type& moment53_rerr ) { + return std::make_pair( isoplanatic_angle(moment53), 3.0/5.0*moment53_rerr ); +} +std::pair< value_type, value_type > time_constant_desi( const value_type& desi, const value_type& desi_rerr ) { + return std::make_pair( time_constant_desi(desi), 3.0/5.0*desi_rerr ); +} +std::pair< value_type, value_type > time_constant( const value_type& wind_moment0, const value_type& wind_moment0_rerr, const value_type& moment0, const value_type& moment0_rerr ) { + return std::make_pair( time_constant( wind_moment0, moment0), std::sqrt(0.1*moment0_rerr*moment0_rerr + 0.5*wind_moment0_rerr*wind_moment0_rerr) ); +} + namespace { +vector_type init_moment_powers() { + vector_type powers(4); + powers(moment_0)=0.0; + powers(moment_53)=5.0/3.0; + powers(moment_1)=1.0; + powers(moment_2)=2.0; + return powers; +} + +vector_type init_wind_moment_powers() { + vector_type powers(1); + powers(moment_0)=0.0; + return powers; +} + // FIXME: vector_type reorder_flat( const atmos::indices::indices_type& indices ) { assert( indices.size1() == atmos::channels ); @@ -51,70 +116,38 @@ } } -const double KM2M = 1000; -const double POWSEE = 0 ; // Power of Cn2 moment for seeing -const double POWISP = 5.0/3.0; // Power of Cn2 moment for isoplanatic angle -const double POWEFF = 1.0; // Power of Cn2 moment for Effective altitude -const double POWISK = 2.0; // Cn2 second moment power - -const double MINRELW = 6e-4 ; // Minimal ratio of the singular value to the maximal one - -static const double powers = { POWSEE, POWISP, POWEFF, POWISK }; -const unsigned int npower = sizeof(powers)/sizeof(powers0); - -static const double wind_powers = { 0.0 }; -const unsigned int wind_npower = sizeof(wind_powers)/sizeof(wind_powers0); - -static const double HBOUND = 0.38;// Height of the boundary layer kilometers to fit the Cn2-integral -static const double HTOP = 25.0; // Height of the hight layer kilometers - -static const int TOPLAYERS = 4; // Number of rejected high-altitude layers -const double JPOW_SEE = 0.6 ; -const double LPOW_SEE = -0.2; -const double LEFF_SEE = 1.1486983549970351; -const double KS = 1.73588e+07; // Calibration constant in seeing formula (AstL 45, p.399): - //\f$ seeing'' = KS * \lambda_{eff}^{-1/5} M_0^{3/5} * - //\sec^{-1}(z)\f$, where \em lambda_eff is in mkm*/ - -const double JPOW_ISP = -0.6; // Power of Cn-5/3th moment for isoplan. angle -const double LPOW_ISP = 1.2; // Power of lambda for isoplan. angle -const double LEFF_ISP = 0.43527528164806206; -const double KP = 0.000756348; // Calibration constant in isoplanatic angle formula (AstL 45, p.399) - // \f$ \theta_0'' = KP * (M_{5/3} /\lambda_{eff}^2)^{-3/5} * \sec^{8/5}(z) \f$, where - // \f$ KP = (2.905*(2\pi)^2*(MKM2M)^{-2})^{-3/5}*206265 \f$*/ -const double KT = 0.175; // Atmospheric time constant empirical calibration constant for DESI for 1,3ms - // integrations in diameter=2cm apertures with -1km shift -const double KL = 0.058 * 2.7464013582652968e-08; - -const double SPOW_TC = -0.6; // Power of DESI for atm. time constant -functor_factory::dcf_matrix_type functor_factory::calculate( const weight_type& wf, unsigned int npower, const double* powers ) { +functor_factory::dcf_matrix_type functor_factory::calculate( const weight_type& wf, const vector_type& powers ) { using namespace boost::numeric::ublas; using namespace utils::algo; using namespace utils::ublas; + const double max_cond = 1667; + const double bound_altitude = 0.38; // Height of the boundary layer kilometers to fit the Cn2-integral + const double top_altitude = 25.0; // Height of the hight layer kilometers + const vector<double> grid( utils::ublas::vector_generator< lg_generator > ( lg_generator(0.0, 0.04, 0.25), 50 ) ); typedef dcf_matrix_type result_type; typedef dcf_matrix_type::size_type size_type; - vector<double>::const_iterator lower = std::find_if( grid.begin(), grid.end(), std::bind2nd( std::greater< double >(), HBOUND) ); - vector<double>::const_reverse_iterator higher = std::find_if( grid.rbegin(), grid.rend(), std::bind2nd( std::less< double >(), HTOP) ); + vector<double>::const_iterator lower = std::find_if( grid.begin(), grid.end(), std::bind2nd( std::greater< double >(), bound_altitude) ); + vector<double>::const_reverse_iterator higher = std::find_if( grid.rbegin(), grid.rend(), std::bind2nd( std::less< double >(), top_altitude) ); /* FIXME: Weights::grid_type::const_iterator */ // both layer altitudes and grid are compile-time contants, so it is assertion. assert( lower != grid.end() ); assert( higher != grid.rend() ); - result_type res( npower, atmos::indices_size ); + result_type res( powers.size(), atmos::indices_size ); vector_range< const vector<double> > h( grid, range( lower.index(), higher.index()+1 ) ); matrix< double > wfv = vector_of_vector_to_matrix( project( grid, wf) ); - for( unsigned int ipow = 0; ipow < npower; ipow++ ) { - double moment_power = powersipow; + for( vector_type::size_type ipow = 0; ipow < powers.size(); ipow++ ) { + double moment_power = powers(ipow); double weight_power = moment_power / 2.0; - double scale = pow( KM2M, moment_power ); + double scale = pow( 1000.0, moment_power ); vector<double> alt = utils::ublas::exp( log( h )*moment_power ) * scale; vector<double> weight = utils::ublas::exp( log( h )*weight_power ); @@ -132,7 +165,7 @@ sd.apply( left, right ); matrix_vector_slice< matrix<double> > singular( diag( A ) ); // it is guaranteed that singular values are sorted - matrix_vector_slice< matrix<double> >::const_iterator small_singular = std::find_if( singular.begin(), singular.end(), std::bind2nd( std::less< double >(), singular(0)*MINRELW) ); + matrix_vector_slice< matrix<double> >::const_iterator small_singular = std::find_if( singular.begin(), singular.end(), std::bind2nd( std::less< double >(), singular(0)/max_cond) ); matrix_vector_slice< matrix<double> >::const_iterator it = singular.begin(); alt = prod( left, alt ); for( ; it != small_singular; ++it ) { @@ -162,16 +195,19 @@ functor_factory::~functor_factory() { } const functor_factory::functor_type& functor_factory::operator() ( const std::string& spectral_type, const weight_type& null_wf, const weight_type& second_wf ) { + static const vector_type moment_powers = init_moment_powers(); + static const vector_type wind_moment_powers = init_wind_moment_powers(); + storage_type::const_iterator it = m_storage.find( spectral_type ); if( it != m_storage.end() ) return it->second; - std::pair< storage_type::iterator, bool > ret = m_storage.insert( std::make_pair( spectral_type, functor_type( calculate( null_wf, npower, powers ), calculate( second_wf, wind_npower, wind_powers ) ) ) ); + std::pair< storage_type::iterator, bool > ret = m_storage.insert( std::make_pair( spectral_type, functor_type( calculate( null_wf, moment_powers ), moment_powers, calculate( second_wf, wind_moment_powers ), wind_moment_powers ) ) ); return ret.first->second; } -functor::result_type functor::operator() ( const indices_type& indices, const indesis_type& indesis, const indices_type& wind_indices, double mz ) const { +functor::result_type functor::operator() ( const indices_type& indices, const indices_type& wind_indices, const indesis_type& desi, double mz ) const { typedef vector< double > vector_type; typedef vector_type::size_type size_type; @@ -183,53 +219,32 @@ aints = prod( m_indices_dcf, aints ); assert( npower == aints.size() ); - for( size_type i = 0; i < npower; ++i ){ - aints(i) = aints(i) / std::pow( mz, 1.0 + powersi );
View file
atmos-2.97.2.tar.gz/src/integrals.h -> atmos-2.97.3.tar.gz/src/integrals.h
Changed
@@ -1,5 +1,5 @@ /* - $Id: integrals.h 200 2011-02-03 14:45:13Z matwey $ + $Id: integrals.h 210 2011-03-02 18:59:06Z matwey $ Copyright (C) 2009 Sternberg Astronomical Institute, MSU This program is free software: you can redistribute it and/or modify @@ -36,26 +36,41 @@ typedef vector< value_type > vector_type; typedef matrix< value_type > matrix_type; + typedef vector_type moments_type; + typedef atmos::indices::indesis_type::value_type desi_type; + + enum moment_power { + moment_0 = 0, + moment_53 = 1, + moment_1 = 2, + moment_2 = 3 + }; + + value_type seeing( const value_type& moment0 ); + value_type effective_altitude( const value_type& moment0, const value_type& moment1 ); + value_type isoplanatic_angle( const value_type& moment53 ); + value_type time_constant_desi( const value_type& desi ); + value_type time_constant( const value_type& wind_moment0, const value_type& moment0 ); + + std::pair< value_type, value_type > seeing( const value_type& moment0, const value_type& moment0_rerr ); + std::pair< value_type, value_type > effective_altitude( const value_type& moment0, const value_type& moment0_rerr, const value_type& moment1, const value_type& moment1_rerr ); + std::pair< value_type, value_type > isoplanatic_angle( const value_type& moment53, const value_type& moment53_rerr ); + std::pair< value_type, value_type > time_constant_desi( const value_type& desi, const value_type& desi_rerr ); + std::pair< value_type, value_type > time_constant( const value_type& wind_moment0, const value_type& wind_moment0_rerr, const value_type& moment0, const value_type& moment0_rerr ); + namespace detail { struct sample { - typedef atmos::integrals::value_type value_type; + typedef atmos::integrals::moments_type moments_type; + typedef atmos::integrals::desi_type desi_type; - value_type free_seeing; - value_type isoplanatic_angle; - value_type altitude_eff; - value_type desi_time_constant; - value_type time_constant; - value_type moment_0; - value_type moment_2; - - sample( const value_type& free_seeing, - const value_type& isoplanatic_angle, - const value_type& altitude_eff, - const value_type& desi_time_constant, - const value_type& time_constant, - const value_type& moment_0, - const value_type& moment_2 ): - free_seeing(free_seeing), isoplanatic_angle(isoplanatic_angle), altitude_eff(altitude_eff), desi_time_constant(desi_time_constant), time_constant(time_constant), moment_0(moment_0), moment_2(moment_2) { + moments_type moments; + moments_type wind_moments; + desi_type desi; + + sample( const moments_type& moments, + const moments_type& wind_moments, + const desi_type& desi ): + moments(moments), wind_moments(wind_moments), desi(desi) { } }; } @@ -64,51 +79,26 @@ public table_storage< detail::sample > { public: typedef storage self_type; -private: - typedef table_storage_column_getter< self_type, detail::sample::value_type > common_column_type; public: - typedef common_column_type free_seeing_column_type; - typedef common_column_type isoplanatic_angle_column_type; - typedef common_column_type altitude_eff_column_type; - typedef common_column_type desi_time_constant_column_type; - typedef common_column_type time_constant_column_type; - typedef common_column_type moment_0_column_type; - typedef common_column_type moment_2_column_type; + typedef table_storage_column_proxy< self_type, detail::sample::moments_type > moments_column_type; + typedef table_storage_column_proxy< self_type, detail::sample::moments_type > wind_moments_column_type; + typedef table_storage_column_proxy< self_type, detail::sample::desi_type > desi_column_type; private: - free_seeing_column_type m_free_seeing_column; - isoplanatic_angle_column_type m_isoplanatic_angle_column; - altitude_eff_column_type m_altitude_eff_column; - desi_time_constant_column_type m_desi_time_constant_column; - time_constant_column_type m_time_constant_column; - moment_0_column_type m_moment_0_column; - moment_2_column_type m_moment_2_column; + moments_column_type m_moments_column; + wind_moments_column_type m_wind_moments_column; + desi_column_type m_desi_column; public: storage( size_type reserve = 60 ); - const free_seeing_column_type& free_seeing() const { return m_free_seeing_column; } - const isoplanatic_angle_column_type& isoplanatic_angle() const { return m_isoplanatic_angle_column; } - const altitude_eff_column_type& altitude_eff() const { return m_altitude_eff_column; } - const desi_time_constant_column_type& desi_time_constant() const { return m_desi_time_constant_column; } - const time_constant_column_type& time_constant() const { return m_time_constant_column; } - const moment_0_column_type& moment_0() const { return m_moment_0_column; } - const moment_2_column_type& moment_2() const { return m_moment_2_column; } + const moments_column_type& moments() const { return m_moments_column; } + const wind_moments_column_type& wind_moments() const { return m_wind_moments_column; } + const desi_column_type& desi() const { return m_desi_column; } - void push_back( const detail::sample& integrals ) { - table_storage< detail::sample >::push_back( integrals ); + void push_back( const moments_type& moments, const moments_type& wind_moments, const desi_type& desi ) { + table_storage< detail::sample >::push_back( detail::sample( moments, wind_moments, desi ) ); } }; -namespace filter { - class regular { - public: - typedef atmos::integrals::detail::sample sample_type; - public: - bool operator() ( const sample_type& x ) const; - private: - bool is_regular( double x ) const; - }; -} - class functor_factory; class functor; @@ -121,7 +111,7 @@ private: storage_type m_storage; - dcf_matrix_type calculate( const weight_type& wf, unsigned int npower, const double* powers ); + dcf_matrix_type calculate( const weight_type& wf, const vector_type& powers ); public: functor_factory(); functor_factory( const functor_factory& factory ); @@ -138,15 +128,18 @@ typedef atmos::indices::indices_type indices_type; typedef atmos::indices::indesis_type indesis_type; - typedef atmos::integrals::detail::sample result_type; + typedef boost::tuple< atmos::integrals::moments_type, atmos::integrals::moments_type, atmos::integrals::desi_type > result_type; private: dcf_matrix_type m_indices_dcf; dcf_matrix_type m_wind_indices_dcf; + + const vector_type& m_powers; + const vector_type& m_wind_powers; public: - functor( const dcf_matrix_type& indices_dcf, const dcf_matrix_type& wind_indices_dcf ): - m_indices_dcf( indices_dcf ), m_wind_indices_dcf( wind_indices_dcf ) {} + functor( const dcf_matrix_type& indices_dcf, const vector_type& powers, const dcf_matrix_type& wind_indices_dcf, const vector_type& wind_powers ): + m_indices_dcf( indices_dcf ), m_powers(powers), m_wind_indices_dcf( wind_indices_dcf ), m_wind_powers(wind_powers) {} - result_type operator() ( const indices_type& indices, const indesis_type& indesis, const indices_type& wind_indices, double mz ) const; + result_type operator() ( const indices_type& indices, const indices_type& wind_indices, const indesis_type& desi, double mz ) const; }; }}
View file
atmos-2.97.2.tar.gz/src/worksheet.cpp -> atmos-2.97.3.tar.gz/src/worksheet.cpp
Changed
@@ -1,5 +1,5 @@ /* - $Id: worksheet.cpp 203 2011-02-03 17:40:12Z matwey $ + $Id: worksheet.cpp 210 2011-03-02 18:59:06Z matwey $ Copyright (C) 2009 Sternberg Astronomical Institute, MSU This program is free software: you can redistribute it and/or modify @@ -219,9 +219,9 @@ double Jtrans_err = dimm::turbulence( context().m_pproxy.dimm_scale(), context().m_pproxy.dimm_Ktransversal(), yerrs(1) ); dimm_turbulence = 0.5*( Jlong + Jtrans ) / m_mean_mz; - dimm_turbulence_err = 0.5 / m_mean_mz * std::sqrt( Jlong_err*Jlong_err + Jtrans_err*Jtrans_err ); + dimm_turbulence_err = std::sqrt( Jlong_err*Jlong_err + Jtrans_err*Jtrans_err ) / (Jlong + Jtrans); dimm_seeing = dimm::seeing( dimm_turbulence ); - dimm_seeing_err = 0.6 * dimm_seeing / dimm_turbulence * dimm_turbulence_err; + dimm_seeing_err = 0.6*dimm_turbulence_err; } indices::flux_type flux_avg = average( m_indices.flux().begin(), m_indices.flux().end(), mean< indices::flux_type >() ); @@ -238,36 +238,39 @@ integrals::storage integrals; for( indices::storage::const_iterator it = m_indices.begin(); it != m_indices.end(); ++it ) { - integrals.push_back( ifun( it->indices, it->indesis, it->wind_indices, m_mean_mz ) ); - } + const boost::tuple<integrals::moments_type,integrals::moments_type,integrals::desi_type>& moments = ifun( it->indices, it->wind_indices, it->indesis, m_mean_mz ); - integrals.unstable_remove_if( integrals::filter::regular() ); + integrals.push_back( moments.get<0>(), moments.get<1>(), moments.get<2>() ); + } - double free_seeing_avg = average( integrals.free_seeing().begin(), integrals.free_seeing().end(), mean< double >() ); - double isoplanatic_angle_avg = average( integrals.isoplanatic_angle().begin(), integrals.isoplanatic_angle().end(), mean< double >() ); - double altitude_eff_avg = average( integrals.altitude_eff().begin(), integrals.altitude_eff().end(), mean< double >() ); - double desi_time_constant_avg= average( integrals.desi_time_constant().begin(),integrals.desi_time_constant().end(),mean< double >() ); - double time_constant_avg = average( integrals.time_constant().begin(), integrals.time_constant().end(), mean< double >() ); - double moment_0_avg = average( integrals.moment_0().begin(), integrals.moment_0().end(), mean< double >() ); - double moment_2_avg = average( integrals.moment_2().begin(), integrals.moment_2().end(), mean< double >() ); - double free_seeing_err = average( integrals.free_seeing().begin(), integrals.free_seeing().end(), error_of_mean< double >(free_seeing_avg) ); - double isoplanatic_angle_err = average( integrals.isoplanatic_angle().begin(), integrals.isoplanatic_angle().end(), error_of_mean< double >(isoplanatic_angle_avg) ); - double altitude_eff_err = average( integrals.altitude_eff().begin(), integrals.altitude_eff().end(), error_of_mean< double >(altitude_eff_avg) ); - double desi_time_constant_err= average( integrals.desi_time_constant().begin(),integrals.desi_time_constant().end(),error_of_mean< double >(desi_time_constant_avg) ); - double time_constant_err = average( integrals.time_constant().begin(), integrals.time_constant().end(), error_of_mean< double >(time_constant_avg) ); - double moment_0_err = average( integrals.moment_0().begin(), integrals.moment_0().end(), error_of_mean< double >(moment_0_avg) ); - double moment_2_err = average( integrals.moment_2().begin(), integrals.moment_2().end(), error_of_mean< double >(moment_2_avg) ); + integrals::moments_type moments = average( integrals.moments().begin(), integrals.moments().end(), mean< integrals::moments_type >() ); + integrals::moments_type moments_err = average( integrals.moments().begin(), integrals.moments().end(), error_of_mean< integrals::moments_type >(moments) ); + integrals::moments_type wind_moments = average( integrals.wind_moments().begin(), integrals.wind_moments().end(), mean< integrals::moments_type >() ); + integrals::moments_type wind_moments_err = average( integrals.wind_moments().begin(), integrals.wind_moments().end(), error_of_mean< integrals::moments_type >(wind_moments) ); + integrals::desi_type desi = average( integrals.desi().begin(), integrals.desi().end(), mean< integrals::desi_type >() ); + integrals::desi_type desi_err = average( integrals.desi().begin(), integrals.desi().end(), error_of_mean< integrals::desi_type >(desi) ); + + double moment_0 = moments(integrals::moment_0); + double moment_2 = moments(integrals::moment_2); + double moment_0_err = moments_err(integrals::moment_0)/moments(integrals::moment_0); + double moment_2_err = moments_err(integrals::moment_2)/moments(integrals::moment_2); + + std::pair<double,double> free_seeing = integrals::seeing(moment_0,moment_0_err); + std::pair<double,double> isoplanatic_angle = integrals::isoplanatic_angle(moments(integrals::moment_53),moments_err(integrals::moment_53)/moments(integrals::moment_53)); + std::pair<double,double> altitude_eff = integrals::effective_altitude(moment_0,moment_0_err,moments(integrals::moment_1),moments_err(integrals::moment_1)/moments(integrals::moment_1)); + std::pair<double,double> desi_time_constant = integrals::time_constant_desi(desi,desi_err/desi); + std::pair<double,double> time_constant = integrals::time_constant(wind_moments(integrals::moment_0),wind_moments_err(integrals::moment_0)/wind_moments(integrals::moment_0),moment_0,moment_0_err); context().m_file_output.atmosphere( meantime, - free_seeing_avg, free_seeing_err, + free_seeing.first, free_seeing.second, dimm_seeing, dimm_seeing_err, dimm_turbulence, dimm_turbulence_err, - isoplanatic_angle_avg, isoplanatic_angle_err, - altitude_eff_avg, altitude_eff_err, - desi_time_constant_avg, desi_time_constant_err, - time_constant_avg, time_constant_err, - moment_0_avg, moment_0_err, - moment_2_avg, moment_2_err ); + isoplanatic_angle.first, isoplanatic_angle.second, + altitude_eff.first, altitude_eff.second, + desi_time_constant.first, desi_time_constant.second, + time_constant.first, time_constant.second, + moment_0, moment_0_err, + moment_2, moment_2_err ); /* mixture */ mixture::storage mixture_samples( m_indices, m_dimm, m_dimm.size() > 1 ? context().m_pproxy.dimm_scale() : 0 );
Locations
Projects
Search
Status Monitor
Help
Open Build Service
OBS Manuals
API Documentation
OBS Portal
Reporting a Bug
Contact
Mailing List
Forums
Chat (IRC)
Twitter
Open Build Service (OBS)
is an
openSUSE project
.