File: model.h

package info (click to toggle)
aoflagger 2.13.0-1
  • links: PTS, VCS
  • area: main
  • in suites: buster
  • size: 4,232 kB
  • sloc: cpp: 61,805; python: 60; sh: 23; makefile: 8
file content (132 lines) | stat: -rw-r--r-- 4,943 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
#ifndef MODEL_H
#define MODEL_H

#include <vector>
#include <cmath>

#include "../structures/image2d.h"

#include "../structures/types.h"

#include "uvimager.h"

template<typename T>
struct OutputReceiver
{
	void SetY(size_t) { }
};
template<>
struct OutputReceiver<UVImager>
{
	UVImager *_imager;
	void SetUVValue(size_t, double u, double v, double r, double i, double w)
	{
		_imager->SetUVValue(u, v, r, i, w);
		_imager->SetUVValue(-u, -v, r, -i, w);
	}
	void SetY(size_t) { }
};
template<>
struct OutputReceiver<TimeFrequencyData>
{
	Image2DPtr _real, _imaginary;
	size_t _y;
	void SetUVValue(size_t x, double, double, double r, double i, double)
	{
		_real->SetValue(x, _y, r);
		_imaginary->SetValue(x, _y, i);
	}
	void SetY(size_t y) { _y = y; }
};

class Model {
	struct Source {
		virtual ~Source() { }
		virtual numl_t Dec(num_t t) const = 0;
		virtual numl_t Ra(num_t t) const = 0;
		virtual numl_t FluxIntensity(num_t t) const = 0;
		virtual numl_t SqrtFluxIntensity(num_t t) const
		{
			return sqrtnl(FluxIntensity(t));
		}
	};
	struct StablePointSource : public Source {
		long double dec, ra, fluxIntensity, sqrtFluxIntensity;
		virtual numl_t Dec(num_t) const final override { return dec; }
		virtual numl_t Ra(num_t) const final override { return ra; }
		virtual numl_t FluxIntensity(num_t) const final override { return fluxIntensity; }
		virtual numl_t SqrtFluxIntensity(num_t) const final override { return sqrtFluxIntensity; }
	};
	struct VariablePointSource : public Source {
		long double dec, ra, fluxIntensity;
		double peakTime, oneOverSigmaSq;
		virtual numl_t Dec(num_t) const final override { return dec; }
		virtual numl_t Ra(num_t) const final override { return ra; }
		virtual numl_t FluxIntensity(num_t t) const final override
		{
			numl_t mu = fmodnl(fabsnl(t-peakTime), 1.0);
			if(mu > 0.5) mu = 1.0 - mu;
			return fluxIntensity * (1.0+expnl(mu*mu*oneOverSigmaSq)) * (1.0 + fmod(t*1007.0, 13.0) / 26.0);
		}
	};
	
	public:
		Model();
		void AddSource(long double dec, long double ra, long double fluxIntensity)
		{
			StablePointSource *source = new StablePointSource();
			source->dec = dec;
			source->ra = ra;
			source->fluxIntensity = fluxIntensity;
			source->sqrtFluxIntensity = sqrt(fluxIntensity);
			_sources.push_back(source);
		}
		void AddVariableSource(long double dec, long double ra, long double fluxIntensity)
		{
			VariablePointSource *source = new VariablePointSource();
			source->dec = dec;
			source->ra = ra;
			source->fluxIntensity = fluxIntensity;
			source->peakTime = 0.2;
			source->oneOverSigmaSq = 1.0/(0.3*0.3);
			_sources.push_back(source);
		}
		void SimulateAntenna(double time, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t frequency, num_t earthLattitude, num_t &r, num_t &i);
		void SimulateUncoherentAntenna(double time, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t frequency, num_t earthLattitude, num_t &r, num_t &i, size_t index);

		template<typename T>
		void SimulateCorrelation(struct OutputReceiver<T> &receiver, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t dz, num_t frequency, num_t channelWidth, double totalTime, double integrationTime);

		void SimulateObservation(class UVImager &imager, class Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA)
		{
			srand(1);
			OutputReceiver<UVImager> imagerOutputter;
			imagerOutputter._imager = &imager;
			SimulateObservation(imagerOutputter, observatorium, delayDirectionDEC, delayDirectionRA);
		}

		std::pair<TimeFrequencyData, TimeFrequencyMetaDataPtr> SimulateObservation(class Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, size_t a1, size_t a2);

		template<typename T>
		void SimulateObservation(struct OutputReceiver<T> &receiver, class Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA);

		static void GetUVPosition(num_t &u, num_t &v, num_t earthLattitudeAngle, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t dz, num_t waveLength);
		static num_t GetWPosition(num_t delayDirectionDec, num_t delayDirectionRA, num_t frequency, num_t earthLattitudeAngle, num_t dx, num_t dy)
		{
			return UVImager::GetWPosition(delayDirectionDec, delayDirectionRA, frequency, earthLattitudeAngle, dx, dy);
		}

		void loadUrsaMajor(double ra, double dec, double factor);
		void loadUrsaMajorDistortingSource(double ra, double dec, double factor, bool slightlyMiss = false);
		void loadUrsaMajorDistortingVariableSource(double ra, double dec, double factor, bool weak=false, bool slightlyMiss=false);
		void loadOnAxisSource(double ra, double dec, double factor);
		
		double NoiseSigma() const { return _noiseSigma; }
		void SetNoiseSigma(double noiseSigma) { _noiseSigma = noiseSigma; }
	private:
		std::vector<Source *> _sources;
		double _noiseSigma, _sourceSigma;
		double _integrationTime;
};

#endif