Struct CloudSystem

Source

pub struct CloudSystem {
    pub layers: Vec<CloudLayer>,
    pub time_offset_s: f64,
}

Fields§

§layers: Vec<CloudLayer>§time_offset_s: f64

Implementations§

Source §

impl CloudSystem

Source

pub fn earth_default() -> Self

Examples found in repository ?

examples/tidal_sim.rs (line 516)

47fn main() {
48    let orbit = EarthOrbit::new();
49    let period_days = orbit.orbital_period_days();
50    let v_perihelion = orbit.velocity_at_distance(orbit.perihelion_m());
51    let v_aphelion = orbit.velocity_at_distance(orbit.aphelion_m());
52    let energy = orbit.specific_orbital_energy();
53    let ang_mom = orbit.specific_angular_momentum();
54    let v_esc = EarthOrbit::escape_velocity_at_surface();
55    let f_sun = orbit.gravitational_force_sun();
56    let r = orbit.current_radius();
57    let v_mean = orbit.mean_orbital_velocity();
58    println!(
59        "Orbit: P={:.2}d  v_peri={:.0}m/s  v_aph={:.0}m/s  v_mean={:.0}m/s",
60        period_days, v_perihelion, v_aphelion, v_mean
61    );
62    println!(
63        "  E={:.2e}J/kg  L={:.2e}m²/s  v_esc={:.0}m/s  F_sun={:.2e}N  r={:.3e}m",
64        energy, ang_mom, v_esc, f_sun, r
65    );
66    println!(
67        "  a={:.3e}m  e={}  i={}°  Ω={}°  ω={}°",
68        SEMI_MAJOR_AXIS,
69        ECCENTRICITY,
70        INCLINATION_DEG,
71        LONGITUDE_ASCENDING_NODE_DEG,
72        ARGUMENT_PERIHELION_DEG
73    );
74
75    let rot = EarthRotation::new();
76    let v_equator = rot.surface_velocity_at_latitude(0.0);
77    let v_45 = rot.surface_velocity_at_latitude(45.0);
78    let accel = rot.centripetal_acceleration_at_latitude(0.0);
79    let coriolis = rot.coriolis_parameter(45.0);
80    let moi = rot.moment_of_inertia();
81    let ke = rot.rotational_kinetic_energy();
82    let l_rot = rot.angular_momentum();
83    let prec = rot.precession_rate_rad_per_year();
84    let daylen_summer = rot.day_length_variation_due_to_tilt(172, 60.0);
85    let daylen_winter = rot.day_length_variation_due_to_tilt(355, 60.0);
86    println!(
87        "Rotation: v_eq={:.0}m/s  v_45={:.0}m/s  a_c={:.4}m/s²  f_45={:.5e}",
88        v_equator, v_45, accel, coriolis
89    );
90    println!(
91        "  I={:.3e}kg·m²  KE={:.3e}J  L={:.3e}kg·m²/s  prec={:.2e}rad/yr",
92        moi, ke, l_rot, prec
93    );
94    println!(
95        "  day@60°N: summer={:.1}h  winter={:.1}h",
96        daylen_summer, daylen_winter
97    );
98    println!(
99        "  sidereal={:.4}s  solar={}s  tilt={}°={:.4}rad  prec_period={}yr",
100        SIDEREAL_DAY_S, SOLAR_DAY_S, AXIAL_TILT_DEG, AXIAL_TILT_RAD, PRECESSION_PERIOD_YEARS
101    );
102
103    let moon_tide = TidalForce::from_moon();
104    let sun_tide = TidalForce::from_sun();
105    let a_moon = moon_tide.tidal_acceleration();
106    let a_sun = sun_tide.tidal_acceleration();
107    let pot = moon_tide.tidal_potential(0.0);
108    let bulge_moon = moon_tide.tidal_bulge_height();
109    let bulge_sun = sun_tide.tidal_bulge_height();
110    let grav = moon_tide.gravitational_attraction();
111    let spring = spring_tide_amplitude();
112    let neap = neap_tide_amplitude();
113    let ratio = lunar_to_solar_tide_ratio();
114    println!(
115        "Tides: a_moon={:.2e}  a_sun={:.2e}  ratio={:.2}",
116        a_moon, a_sun, ratio
117    );
118    println!(
119        "  bulge_moon={:.3}m  bulge_sun={:.3}m  spring={:.3}m  neap={:.3}m",
120        bulge_moon, bulge_sun, spring, neap
121    );
122    println!("  potential(0)={:.2e}  F_grav={:.2e}N", pot, grav);
123
124    let chix = chicxulub_equivalent();
125    let tung = tunguska_equivalent();
126    let ke_chix = chix.kinetic_energy_mt();
127    let crater = chix.crater_diameter_m(2700.0);
128    let fireball = chix.fireball_radius_m();
129    let ejecta = chix.ejecta_volume_m3(2700.0);
130    let v_imp = chix.impact_velocity();
131    let ke_tung = tung.kinetic_energy_mt();
132    println!(
133        "Chicxulub: {:.2e}Mt  crater={:.0}m  fireball={:.0}m  ejecta={:.2e}m³  v_imp={:.0}m/s",
134        ke_chix, crater, fireball, ejecta, v_imp
135    );
136    println!("Tunguska: {:.2e}Mt", ke_tung);
137
138    let mut moon = MoonState::new();
139    let source = match moon.source.get() {
140        MoonSource::Binary => "binary",
141        MoonSource::Simulation => "simulated",
142    };
143    let pos0 = moon.position();
144    moon.step(3600.0);
145    let pos1 = moon.position();
146    let g_moon = moon.gravity_at(EARTH_MOON_DISTANCE);
147    println!(
148        "Moon({}): pos0=({:.0},{:.0})  pos1h=({:.0},{:.0})  g@dist={:.4e}m/s²",
149        source, pos0.0, pos0.1, pos1.0, pos1.1, g_moon
150    );
151    println!(
152        "  M={:.3e}kg  R={:.4e}m  d={:.3e}m  P={:.0}s",
153        LUNAR_MASS, LUNAR_RADIUS, EARTH_MOON_DISTANCE, LUNAR_ORBITAL_PERIOD
154    );
155
156    let iss = ArtificialSatellite::leo("ISS", 420_000.0, 408_000.0);
157    let geo = ArtificialSatellite::geo("GEO-Sat", 300.0);
158    let custom = ArtificialSatellite::new("Molniya", 1500.0, 500_000.0, 0.7, 1.1);
159    let p_iss = iss.orbital_period_s();
160    let v_iss = iss.orbital_velocity_ms();
161    let p_geo = geo.orbital_period_s();
162    let g_surf = iss.gravity_at_surface();
163    let pos_iss = iss.position();
164    let pos_geo = geo.position();
165
166    let orbit_type = |sat: &ArtificialSatellite| -> &str {
167        match sat.orbit_type {
168            OrbitType::LEO => "LEO",
169            OrbitType::MEO => "MEO",
170            OrbitType::GEO => "GEO",
171            OrbitType::HEO => "HEO",
172            OrbitType::Custom => "Custom",
173        }
174    };
175
176    println!(
177        "ISS({}): P={:.0}s  v={:.0}m/s  g_surf={:.2}  pos=({:.0},{:.0},{:.0})",
178        orbit_type(&iss),
179        p_iss,
180        v_iss,
181        g_surf,
182        pos_iss.0,
183        pos_iss.1,
184        pos_iss.2
185    );
186    println!(
187        "GEO({}): P={:.0}s  pos=({:.0},{:.0},{:.0})",
188        orbit_type(&geo),
189        p_geo,
190        pos_geo.0,
191        pos_geo.1,
192        pos_geo.2
193    );
194    println!(
195        "Custom({}): e={:.1}  i={:.1}rad",
196        orbit_type(&custom),
197        custom.eccentricity,
198        custom.inclination_rad
199    );
200
201    let mut constellation = Constellation::new("Starlink-test");
202    constellation.add(ArtificialSatellite::leo("S1", 260.0, 550_000.0));
203    constellation.add(ArtificialSatellite::leo("S2", 260.0, 550_000.0));
204    let mut sat_stepped = iss;
205    sat_stepped.step(60.0);
206    constellation.step_all(60.0);
207    let positions = constellation.positions();
208    println!(
209        "Constellation: {} sats  positions={}",
210        positions.len(),
211        positions.len()
212    );
213
214    let dt = DateTime::new(2024, 6, 21, 12, 0, 0.0);
215    let jd = dt.to_julian_date();
216    let dt_back = DateTime::from_julian_date(jd);
217    let unix = dt.to_unix_timestamp();
218    let dt_unix = DateTime::from_unix_timestamp(unix);
219    println!(
220        "DateTime: {}-{}-{} {}:{}:{:.0}  JD={:.4}  unix={:.0}",
221        dt.year, dt.month, dt.day, dt.hour, dt.minute, dt.second, jd, unix
222    );
223    println!(
224        "  roundtrip JD: {}-{}-{}  roundtrip unix: {}-{}-{}",
225        dt_back.year, dt_back.month, dt_back.day, dt_unix.year, dt_unix.month, dt_unix.day
226    );
227    println!(
228        "  J2000_JD={} UNIX_JD={} SPD={}",
229        J2000_JD, UNIX_EPOCH_JD, SECONDS_PER_DAY
230    );
231
232    let mut epoch = Epoch::j2000();
233    let epoch_mjd = Epoch::from_mjd(51544.5);
234    let centuries = epoch.centuries_since_j2000();
235    let days = epoch.days_since_j2000();
236    let gmst = epoch.gmst_degrees();
237    let mjd = epoch.to_mjd();
238    epoch.advance_days(365.25);
239    let centuries_1yr = epoch.centuries_since_j2000();
240    let mut epoch2 = Epoch::from_jd(J2000_JD);
241    epoch2.advance_seconds(86400.0);
242    println!(
243        "Epoch: T0={:.4}  days={:.1}  GMST={:.2}°  MJD={:.1}",
244        centuries, days, gmst, mjd
245    );
246    println!(
247        "  +1yr: T={:.6}  from_mjd.jd={:.1}  epoch2.days={:.1}",
248        centuries_1yr,
249        epoch_mjd.julian_date,
250        epoch2.days_since_j2000()
251    );
252    println!(
253        "  J2000={} J1950={} MJD_OFF={}",
254        J2000_EPOCH, J1950_EPOCH, MJD_OFFSET
255    );
256
257    let mut ts = TimeScale::realtime();
258    ts.step(1.0);
259    let sim_dt = ts.simulation_dt(1.0);
260    let hours = ts.sim_hours();
261    let days_ts = ts.sim_days();
262    let years = ts.sim_years();
263    ts.pause();
264    ts.resume();
265    ts.toggle_pause();
266    ts.set_speed(100.0);
267    let mut ff = TimeScale::fast_forward(10.0);
268    ff.step(1.0);
269    let sm = TimeScale::slow_motion(2.0);
270    println!(
271        "TimeScale: dt={:.1}  h={:.6}  d={:.8}  yr={:.10}  ff_sim={:.1}s  sm_speed={:.1}",
272        sim_dt, hours, days_ts, years, ff.simulation_time_s, sm.speed_multiplier
273    );
274
275    let sun_pos = SolarPosition::compute(jd, 48.8, 2.3);
276    let above = sun_pos.is_above_horizon();
277    let dist_sun = sun_pos.distance_m();
278    println!(
279        "Solar: el={:.1}°  az={:.1}°  above={}  dist={:.3e}m  tilt={}°",
280        sun_pos.elevation_deg, sun_pos.azimuth_deg, above, dist_sun, EARTH_AXIAL_TILT_DEG
281    );
282
283    let dnc = DayNightCycle::new(jd);
284    let state_paris = dnc.state_at(48.8, 2.3);
285    let state_text = match state_paris {
286        DaylightState::Day => "day",
287        DaylightState::Night => "night",
288        DaylightState::CivilTwilight => "civil_twilight",
289        DaylightState::NauticalTwilight => "nautical_twilight",
290        DaylightState::AstronomicalTwilight => "astro_twilight",
291    };
292    let terminator = dnc.terminator_points(36);
293    let ambient = dnc.ambient_light(48.8, 2.3);
294    println!(
295        "DayNight: Paris={}  ambient={:.2}  terminator_pts={}",
296        state_text,
297        ambient,
298        terminator.len()
299    );
300
301    let season_state = season_at(jd, 48.8);
302    let season_name = match season_state.season_north {
303        Season::Spring => "spring",
304        Season::Summer => "summer",
305        Season::Autumn => "autumn",
306        Season::Winter => "winter",
307    };
308    let south_name = match season_state.season_south {
309        Season::Spring => "spring",
310        Season::Summer => "summer",
311        Season::Autumn => "autumn",
312        Season::Winter => "winter",
313    };
314    let (sub_lat, sub_lon) = subsolar_point(jd);
315    println!(
316        "Season: N={}  S={}  decl={:.2}°  day_h={:.2}  subsolar=({:.1},{:.1})",
317        season_name,
318        south_name,
319        season_state.solar_declination_deg,
320        season_state.day_length_hours,
321        sub_lat,
322        sub_lon
323    );
324    println!(
325        "  tilt={}°  trop_year={}d  vernal_jd={}",
326        SEASONS_TILT, TROPICAL_YEAR_DAYS, VERNAL_EQUINOX_JD
327    );
328
329    let paris = LatLon::new(48.8566, 2.3522, 35.0);
330    let tokyo = LatLon::new(35.6762, 139.6503, 40.0);
331    let ecef_paris = paris.to_ecef();
332    let cart = paris.to_cartesian_simple();
333    let dist_pt = paris.distance_to(&tokyo);
334    let latlon_back = ecef_paris.to_latlon();
335    println!(
336        "Paris->ECEF: ({:.0},{:.0},{:.0})  cart=({:.0},{:.0},{:.0})",
337        ecef_paris.x, ecef_paris.y, ecef_paris.z, cart[0], cart[1], cart[2]
338    );
339    println!(
340        "  dist(Paris-Tokyo)={:.0}m  roundtrip_lat={:.4}°",
341        dist_pt, latlon_back.lat_deg
342    );
343    println!(
344        "  f={:.9}  a={:.0}m  b={:.3}m",
345        EARTH_FLATTENING, EARTH_SEMI_MAJOR_M, EARTH_SEMI_MINOR_M
346    );
347
348    let regions = RegionDatabase::continents();
349    let region = regions.point_in_region(48.8, 2.3);
350    if let Some(r) = region {
351        let rtype = match r.region_type {
352            RegionType::Continent => "continent",
353            RegionType::Ocean => "ocean",
354            RegionType::Sea => "sea",
355            RegionType::Country => "country",
356            RegionType::Island => "island",
357        };
358        println!(
359            "Region(48.8,2.3): {} ({})  area={:.0}km²",
360            r.name, rtype, r.area_km2
361        );
362    }
363
364    let elev_prov = ElevationProvider::global(30.0);
365    let elev_everest = elev_prov.sample(27.988, 86.925);
366    let notables = ElevationProvider::notable_elevations();
367    println!(
368        "Elevation(Everest)={:.0}m  notables={}",
369        elev_everest,
370        notables.len()
371    );
372
373    let bathy = BathymetryData::global(60.0);
374    let depth_mariana = bathy.sample(11.35, 142.2);
375    let is_ocean = bathy.is_ocean(0.0, -30.0);
376    let stats = BathymetryData::ocean_stats();
377    let basins = OceanBasin::major_basins();
378    println!(
379        "Bathy(Mariana)={:.0}m  is_ocean(0,-30)={}  basins={}  avg_depth={:.0}m",
380        depth_mariana,
381        is_ocean,
382        basins.len(),
383        stats.avg_depth_m
384    );
385
386    let mut hmap = Heightmap::new(64);
387    hmap.generate_procedural(42, 6, 0.5, 2.0);
388    let h_sample = hmap.sample(45.0, 10.0);
389    let r_at = hmap.radius_at(45.0, 10.0);
390    println!(
391        "Heightmap(64): sample(45,10)={:.1}m  radius={:.0}m",
392        h_sample, r_at
393    );
394
395    let config = LodConfig::default();
396    let mut lod = LodTerrain::new(config);
397    lod.update([6_371_000.0 + 1000.0, 0.0, 0.0]);
398
399    let mesh = TerrainMesh::from_region(0.0, 10.0, 0.0, 10.0, 8, &|lat, lon| hmap.sample(lat, lon));
400    println!(
401        "Mesh: {} vertices  {} triangles",
402        mesh.vertex_count(),
403        mesh.triangle_count()
404    );
405
406    let classifier = BiomeClassifier::default();
407    let biome = classifier.classify(2000.0, 45.0, 0.6);
408    let biome_name = match biome {
409        Biome::Ocean => "ocean",
410        Biome::Beach => "beach",
411        Biome::Desert => "desert",
412        Biome::Grassland => "grassland",
413        Biome::Forest => "forest",
414        Biome::Tundra => "tundra",
415        Biome::Snow => "snow",
416        Biome::Mountain => "mountain",
417        Biome::Volcanic => "volcanic",
418        Biome::Taiga => "taiga",
419    };
420    let splat = classifier.splat(2000.0, 45.0, 0.6);
421    println!(
422        "Biome(2000m,45°,0.6)={}  splat_top={:.2}",
423        biome_name, splat.weights[0].1
424    );
425
426    let face = Face::PosZ;
427    println!("LOD face={:?}", face);
428
429    let mat_ocean = PbrMaterial::ocean();
430    let mat_grass = PbrMaterial::grassland();
431    let mat_desert = PbrMaterial::desert();
432    let mat_snow = PbrMaterial::snow();
433    let mat_rock = PbrMaterial::rock();
434    let mat_volc = PbrMaterial::volcanic();
435    let mat_forest = PbrMaterial::forest();
436    let mat_ice = PbrMaterial::ice();
437    println!(
438        "Materials: ocean_rough={:.2}  grass_rough={:.2}  desert_metal={:.2}  \
439        snow_rough={:.2}  rock={:.2}  volc={:.2}  forest={:.2}  ice={:.2}",
440        mat_ocean.roughness,
441        mat_grass.roughness,
442        mat_desert.metallic,
443        mat_snow.roughness,
444        mat_rock.roughness,
445        mat_volc.roughness,
446        mat_forest.roughness,
447        mat_ice.roughness
448    );
449
450    let sh_terrain = ShaderData::terrain();
451    let sh_atm = ShaderData::atmosphere();
452    let sh_ocean = ShaderData::ocean();
453    println!(
454        "Shaders: terrain={} uniforms  atm={}  ocean={}",
455        sh_terrain.uniforms.len(),
456        sh_atm.uniforms.len(),
457        sh_ocean.uniforms.len()
458    );
459    for u in &sh_terrain.uniforms {
460        let val = match &u.value {
461            UniformValue::Float(f) => format!("f{:.2}", f),
462            UniformValue::Vec3(v) => format!("v3({:.1},{:.1},{:.1})", v[0], v[1], v[2]),
463            UniformValue::Vec4(v) => format!("v4({:.1},{:.1},{:.1},{:.1})", v[0], v[1], v[2], v[3]),
464            UniformValue::Mat4(_) => "mat4".to_string(),
465            UniformValue::Int(i) => format!("i{}", i),
466        };
467        println!("  {}={}", u.name, val);
468    }
469
470    let atm_params = AtmosphereParams::default();
471    let ray_dens = atm_params.rayleigh_density(4000.0);
472    let mie_dens = atm_params.mie_density(600.0);
473    let scat_r = atm_params.scatter_rayleigh(0.5);
474    let scat_m = atm_params.scatter_mie(0.5);
475    let sky = atm_params.sky_color([0.0, 1.0, 0.0], [0.0, 0.5, 0.866], 0.0);
476    println!(
477        "Atmosphere: ρ_ray(4km)={:.4}  ρ_mie(600m)={:.4}  scat_r={:.4}  scat_m={:.4}",
478        ray_dens, mie_dens, scat_r, scat_m
479    );
480    println!("  sky=({:.4},{:.4},{:.4})", sky[0], sky[1], sky[2]);
481
482    let ocean_rdr = OceanParams::default();
483    let phillips = ocean_rdr.phillips_spectrum(0.1, 0.05);
484    let spectrum = ocean_rdr.generate_spectrum();
485    let disp = ocean_rdr.dispersion(0.1);
486    println!(
487        "Ocean render: phillips(0.1,0.05)={:.6}  spectrum_size={}  disp={:.4}",
488        phillips,
489        spectrum.heights.len(),
490        disp
491    );
492
493    let cumulus = CloudLayer::cumulus();
494    let stratus = CloudLayer::stratus();
495    let cirrus = CloudLayer::cirrus();
496    let cb = CloudLayer::cumulonimbus();
497    let cloud_type = |ct: &CloudType| match ct {
498        CloudType::Cumulus => "Cu",
499        CloudType::Stratus => "St",
500        CloudType::Cirrus => "Ci",
501        CloudType::Cumulonimbus => "Cb",
502        CloudType::Altocumulus => "Ac",
503    };
504    println!(
505        "Clouds: {} base={:.0}m  {} base={:.0}m  {} base={:.0}m  {} base={:.0}m",
506        cloud_type(&cumulus.cloud_type),
507        cumulus.base_altitude_m,
508        cloud_type(&stratus.cloud_type),
509        stratus.base_altitude_m,
510        cloud_type(&cirrus.cloud_type),
511        cirrus.base_altitude_m,
512        cloud_type(&cb.cloud_type),
513        cb.base_altitude_m
514    );
515
516    let mut clouds = CloudSystem::earth_default();
517    clouds.step(3600.0);
518    let density = clouds.sample_density(5000.0, 45.0, 10.0);
519    println!("CloudSystem: density@5km={:.4}", density);
520
521    let rainforest = tropical_rainforest();
522    let boreal = boreal_forest();
523    let reef = coral_reef();
524    println!("Ecosystems:");
525    for eco in [&rainforest, &boreal, &reef] {
526        let sp = eco.total_species();
527        let ind = eco.total_individuals();
528        let shannon = eco.shannon_index();
529        let simpson = eco.simpson_diversity();
530        let area_sp = eco.expected_species_from_area(0.25, 100.0);
531        let npp = eco.total_npp_gc_yr();
532        let turnover = eco.biomass_turnover_time_yr(15_000.0);
533        println!(
534            "  {} species  {} ind  H'={:.2}  D={:.4}  SAR={:.0}  NPP={:.0}  τ={:.1}yr",
535            sp, ind, shannon, simpson, area_sp, npp, turnover
536        );
537    }
538
539    let tb = tropical_broadleaf();
540    let tg = temperate_grassland();
541    let cf = coniferous_forest();
542    let photo = tb.photosynthesis_rate(400.0, 1500.0, 28.0);
543    let canopy = tb.canopy_photosynthesis(photo);
544    let transp = tg.transpiration_mm_day(1.5, 0.3);
545    let npp_veg = cf.npp_kgc_m2_yr(15.0);
546    let crt = cf.carbon_residence_time_yr(npp_veg);
547    let light = beer_lambert_light_extinction(2000.0, 4.0, 0.5);
548    println!(
549        "Vegetation: photo={:.2}  canopy={:.2}  transp={:.2}mm/d  NPP={:.3}  CRT={:.1}yr  light={:.1}",
550        photo, canopy, transp, npp_veg, crt, light
551    );
552
553    let elephant = african_elephant();
554    let whale = blue_whale();
555    let gr_el = elephant.growth_rate();
556    let proj = elephant.project_forward(10.0);
557    let metab = whale.metabolic_rate_w();
558    let range = elephant.home_range_km2();
559    let gen_time = whale.generation_time_years();
560    let lifespan = elephant.max_lifespan_years();
561    println!(
562        "Elephant: r={:.3}  N(+10yr)={:.0}  range={:.1}km²  lifespan={:.0}yr",
563        gr_el, proj, range, lifespan
564    );
565    println!("Whale: metab={:.0}W  gen={:.1}yr", metab, gen_time);
566
567    let mut pp = PredatorPrey {
568        prey: african_elephant(),
569        predator: blue_whale(),
570        attack_rate: 0.01,
571        conversion_efficiency: 0.1,
572        predator_death_rate: 0.05,
573    };
574    let prey_r = pp.prey_growth_rate();
575    let pred_r = pp.predator_growth_rate();
576    pp.step(0.1);
577    println!(
578        "Lotka-Volterra: prey_r={:.1}  pred_r={:.1}  -> prey={:.1}  pred={:.2}",
579        prey_r, pred_r, pp.prey.count, pp.predator.count
580    );
581}

Source

pub fn step(&mut self, dt_s: f64)

Examples found in repository ?

examples/tidal_sim.rs (line 517)

47fn main() {
48    let orbit = EarthOrbit::new();
49    let period_days = orbit.orbital_period_days();
50    let v_perihelion = orbit.velocity_at_distance(orbit.perihelion_m());
51    let v_aphelion = orbit.velocity_at_distance(orbit.aphelion_m());
52    let energy = orbit.specific_orbital_energy();
53    let ang_mom = orbit.specific_angular_momentum();
54    let v_esc = EarthOrbit::escape_velocity_at_surface();
55    let f_sun = orbit.gravitational_force_sun();
56    let r = orbit.current_radius();
57    let v_mean = orbit.mean_orbital_velocity();
58    println!(
59        "Orbit: P={:.2}d  v_peri={:.0}m/s  v_aph={:.0}m/s  v_mean={:.0}m/s",
60        period_days, v_perihelion, v_aphelion, v_mean
61    );
62    println!(
63        "  E={:.2e}J/kg  L={:.2e}m²/s  v_esc={:.0}m/s  F_sun={:.2e}N  r={:.3e}m",
64        energy, ang_mom, v_esc, f_sun, r
65    );
66    println!(
67        "  a={:.3e}m  e={}  i={}°  Ω={}°  ω={}°",
68        SEMI_MAJOR_AXIS,
69        ECCENTRICITY,
70        INCLINATION_DEG,
71        LONGITUDE_ASCENDING_NODE_DEG,
72        ARGUMENT_PERIHELION_DEG
73    );
74
75    let rot = EarthRotation::new();
76    let v_equator = rot.surface_velocity_at_latitude(0.0);
77    let v_45 = rot.surface_velocity_at_latitude(45.0);
78    let accel = rot.centripetal_acceleration_at_latitude(0.0);
79    let coriolis = rot.coriolis_parameter(45.0);
80    let moi = rot.moment_of_inertia();
81    let ke = rot.rotational_kinetic_energy();
82    let l_rot = rot.angular_momentum();
83    let prec = rot.precession_rate_rad_per_year();
84    let daylen_summer = rot.day_length_variation_due_to_tilt(172, 60.0);
85    let daylen_winter = rot.day_length_variation_due_to_tilt(355, 60.0);
86    println!(
87        "Rotation: v_eq={:.0}m/s  v_45={:.0}m/s  a_c={:.4}m/s²  f_45={:.5e}",
88        v_equator, v_45, accel, coriolis
89    );
90    println!(
91        "  I={:.3e}kg·m²  KE={:.3e}J  L={:.3e}kg·m²/s  prec={:.2e}rad/yr",
92        moi, ke, l_rot, prec
93    );
94    println!(
95        "  day@60°N: summer={:.1}h  winter={:.1}h",
96        daylen_summer, daylen_winter
97    );
98    println!(
99        "  sidereal={:.4}s  solar={}s  tilt={}°={:.4}rad  prec_period={}yr",
100        SIDEREAL_DAY_S, SOLAR_DAY_S, AXIAL_TILT_DEG, AXIAL_TILT_RAD, PRECESSION_PERIOD_YEARS
101    );
102
103    let moon_tide = TidalForce::from_moon();
104    let sun_tide = TidalForce::from_sun();
105    let a_moon = moon_tide.tidal_acceleration();
106    let a_sun = sun_tide.tidal_acceleration();
107    let pot = moon_tide.tidal_potential(0.0);
108    let bulge_moon = moon_tide.tidal_bulge_height();
109    let bulge_sun = sun_tide.tidal_bulge_height();
110    let grav = moon_tide.gravitational_attraction();
111    let spring = spring_tide_amplitude();
112    let neap = neap_tide_amplitude();
113    let ratio = lunar_to_solar_tide_ratio();
114    println!(
115        "Tides: a_moon={:.2e}  a_sun={:.2e}  ratio={:.2}",
116        a_moon, a_sun, ratio
117    );
118    println!(
119        "  bulge_moon={:.3}m  bulge_sun={:.3}m  spring={:.3}m  neap={:.3}m",
120        bulge_moon, bulge_sun, spring, neap
121    );
122    println!("  potential(0)={:.2e}  F_grav={:.2e}N", pot, grav);
123
124    let chix = chicxulub_equivalent();
125    let tung = tunguska_equivalent();
126    let ke_chix = chix.kinetic_energy_mt();
127    let crater = chix.crater_diameter_m(2700.0);
128    let fireball = chix.fireball_radius_m();
129    let ejecta = chix.ejecta_volume_m3(2700.0);
130    let v_imp = chix.impact_velocity();
131    let ke_tung = tung.kinetic_energy_mt();
132    println!(
133        "Chicxulub: {:.2e}Mt  crater={:.0}m  fireball={:.0}m  ejecta={:.2e}m³  v_imp={:.0}m/s",
134        ke_chix, crater, fireball, ejecta, v_imp
135    );
136    println!("Tunguska: {:.2e}Mt", ke_tung);
137
138    let mut moon = MoonState::new();
139    let source = match moon.source.get() {
140        MoonSource::Binary => "binary",
141        MoonSource::Simulation => "simulated",
142    };
143    let pos0 = moon.position();
144    moon.step(3600.0);
145    let pos1 = moon.position();
146    let g_moon = moon.gravity_at(EARTH_MOON_DISTANCE);
147    println!(
148        "Moon({}): pos0=({:.0},{:.0})  pos1h=({:.0},{:.0})  g@dist={:.4e}m/s²",
149        source, pos0.0, pos0.1, pos1.0, pos1.1, g_moon
150    );
151    println!(
152        "  M={:.3e}kg  R={:.4e}m  d={:.3e}m  P={:.0}s",
153        LUNAR_MASS, LUNAR_RADIUS, EARTH_MOON_DISTANCE, LUNAR_ORBITAL_PERIOD
154    );
155
156    let iss = ArtificialSatellite::leo("ISS", 420_000.0, 408_000.0);
157    let geo = ArtificialSatellite::geo("GEO-Sat", 300.0);
158    let custom = ArtificialSatellite::new("Molniya", 1500.0, 500_000.0, 0.7, 1.1);
159    let p_iss = iss.orbital_period_s();
160    let v_iss = iss.orbital_velocity_ms();
161    let p_geo = geo.orbital_period_s();
162    let g_surf = iss.gravity_at_surface();
163    let pos_iss = iss.position();
164    let pos_geo = geo.position();
165
166    let orbit_type = |sat: &ArtificialSatellite| -> &str {
167        match sat.orbit_type {
168            OrbitType::LEO => "LEO",
169            OrbitType::MEO => "MEO",
170            OrbitType::GEO => "GEO",
171            OrbitType::HEO => "HEO",
172            OrbitType::Custom => "Custom",
173        }
174    };
175
176    println!(
177        "ISS({}): P={:.0}s  v={:.0}m/s  g_surf={:.2}  pos=({:.0},{:.0},{:.0})",
178        orbit_type(&iss),
179        p_iss,
180        v_iss,
181        g_surf,
182        pos_iss.0,
183        pos_iss.1,
184        pos_iss.2
185    );
186    println!(
187        "GEO({}): P={:.0}s  pos=({:.0},{:.0},{:.0})",
188        orbit_type(&geo),
189        p_geo,
190        pos_geo.0,
191        pos_geo.1,
192        pos_geo.2
193    );
194    println!(
195        "Custom({}): e={:.1}  i={:.1}rad",
196        orbit_type(&custom),
197        custom.eccentricity,
198        custom.inclination_rad
199    );
200
201    let mut constellation = Constellation::new("Starlink-test");
202    constellation.add(ArtificialSatellite::leo("S1", 260.0, 550_000.0));
203    constellation.add(ArtificialSatellite::leo("S2", 260.0, 550_000.0));
204    let mut sat_stepped = iss;
205    sat_stepped.step(60.0);
206    constellation.step_all(60.0);
207    let positions = constellation.positions();
208    println!(
209        "Constellation: {} sats  positions={}",
210        positions.len(),
211        positions.len()
212    );
213
214    let dt = DateTime::new(2024, 6, 21, 12, 0, 0.0);
215    let jd = dt.to_julian_date();
216    let dt_back = DateTime::from_julian_date(jd);
217    let unix = dt.to_unix_timestamp();
218    let dt_unix = DateTime::from_unix_timestamp(unix);
219    println!(
220        "DateTime: {}-{}-{} {}:{}:{:.0}  JD={:.4}  unix={:.0}",
221        dt.year, dt.month, dt.day, dt.hour, dt.minute, dt.second, jd, unix
222    );
223    println!(
224        "  roundtrip JD: {}-{}-{}  roundtrip unix: {}-{}-{}",
225        dt_back.year, dt_back.month, dt_back.day, dt_unix.year, dt_unix.month, dt_unix.day
226    );
227    println!(
228        "  J2000_JD={} UNIX_JD={} SPD={}",
229        J2000_JD, UNIX_EPOCH_JD, SECONDS_PER_DAY
230    );
231
232    let mut epoch = Epoch::j2000();
233    let epoch_mjd = Epoch::from_mjd(51544.5);
234    let centuries = epoch.centuries_since_j2000();
235    let days = epoch.days_since_j2000();
236    let gmst = epoch.gmst_degrees();
237    let mjd = epoch.to_mjd();
238    epoch.advance_days(365.25);
239    let centuries_1yr = epoch.centuries_since_j2000();
240    let mut epoch2 = Epoch::from_jd(J2000_JD);
241    epoch2.advance_seconds(86400.0);
242    println!(
243        "Epoch: T0={:.4}  days={:.1}  GMST={:.2}°  MJD={:.1}",
244        centuries, days, gmst, mjd
245    );
246    println!(
247        "  +1yr: T={:.6}  from_mjd.jd={:.1}  epoch2.days={:.1}",
248        centuries_1yr,
249        epoch_mjd.julian_date,
250        epoch2.days_since_j2000()
251    );
252    println!(
253        "  J2000={} J1950={} MJD_OFF={}",
254        J2000_EPOCH, J1950_EPOCH, MJD_OFFSET
255    );
256
257    let mut ts = TimeScale::realtime();
258    ts.step(1.0);
259    let sim_dt = ts.simulation_dt(1.0);
260    let hours = ts.sim_hours();
261    let days_ts = ts.sim_days();
262    let years = ts.sim_years();
263    ts.pause();
264    ts.resume();
265    ts.toggle_pause();
266    ts.set_speed(100.0);
267    let mut ff = TimeScale::fast_forward(10.0);
268    ff.step(1.0);
269    let sm = TimeScale::slow_motion(2.0);
270    println!(
271        "TimeScale: dt={:.1}  h={:.6}  d={:.8}  yr={:.10}  ff_sim={:.1}s  sm_speed={:.1}",
272        sim_dt, hours, days_ts, years, ff.simulation_time_s, sm.speed_multiplier
273    );
274
275    let sun_pos = SolarPosition::compute(jd, 48.8, 2.3);
276    let above = sun_pos.is_above_horizon();
277    let dist_sun = sun_pos.distance_m();
278    println!(
279        "Solar: el={:.1}°  az={:.1}°  above={}  dist={:.3e}m  tilt={}°",
280        sun_pos.elevation_deg, sun_pos.azimuth_deg, above, dist_sun, EARTH_AXIAL_TILT_DEG
281    );
282
283    let dnc = DayNightCycle::new(jd);
284    let state_paris = dnc.state_at(48.8, 2.3);
285    let state_text = match state_paris {
286        DaylightState::Day => "day",
287        DaylightState::Night => "night",
288        DaylightState::CivilTwilight => "civil_twilight",
289        DaylightState::NauticalTwilight => "nautical_twilight",
290        DaylightState::AstronomicalTwilight => "astro_twilight",
291    };
292    let terminator = dnc.terminator_points(36);
293    let ambient = dnc.ambient_light(48.8, 2.3);
294    println!(
295        "DayNight: Paris={}  ambient={:.2}  terminator_pts={}",
296        state_text,
297        ambient,
298        terminator.len()
299    );
300
301    let season_state = season_at(jd, 48.8);
302    let season_name = match season_state.season_north {
303        Season::Spring => "spring",
304        Season::Summer => "summer",
305        Season::Autumn => "autumn",
306        Season::Winter => "winter",
307    };
308    let south_name = match season_state.season_south {
309        Season::Spring => "spring",
310        Season::Summer => "summer",
311        Season::Autumn => "autumn",
312        Season::Winter => "winter",
313    };
314    let (sub_lat, sub_lon) = subsolar_point(jd);
315    println!(
316        "Season: N={}  S={}  decl={:.2}°  day_h={:.2}  subsolar=({:.1},{:.1})",
317        season_name,
318        south_name,
319        season_state.solar_declination_deg,
320        season_state.day_length_hours,
321        sub_lat,
322        sub_lon
323    );
324    println!(
325        "  tilt={}°  trop_year={}d  vernal_jd={}",
326        SEASONS_TILT, TROPICAL_YEAR_DAYS, VERNAL_EQUINOX_JD
327    );
328
329    let paris = LatLon::new(48.8566, 2.3522, 35.0);
330    let tokyo = LatLon::new(35.6762, 139.6503, 40.0);
331    let ecef_paris = paris.to_ecef();
332    let cart = paris.to_cartesian_simple();
333    let dist_pt = paris.distance_to(&tokyo);
334    let latlon_back = ecef_paris.to_latlon();
335    println!(
336        "Paris->ECEF: ({:.0},{:.0},{:.0})  cart=({:.0},{:.0},{:.0})",
337        ecef_paris.x, ecef_paris.y, ecef_paris.z, cart[0], cart[1], cart[2]
338    );
339    println!(
340        "  dist(Paris-Tokyo)={:.0}m  roundtrip_lat={:.4}°",
341        dist_pt, latlon_back.lat_deg
342    );
343    println!(
344        "  f={:.9}  a={:.0}m  b={:.3}m",
345        EARTH_FLATTENING, EARTH_SEMI_MAJOR_M, EARTH_SEMI_MINOR_M
346    );
347
348    let regions = RegionDatabase::continents();
349    let region = regions.point_in_region(48.8, 2.3);
350    if let Some(r) = region {
351        let rtype = match r.region_type {
352            RegionType::Continent => "continent",
353            RegionType::Ocean => "ocean",
354            RegionType::Sea => "sea",
355            RegionType::Country => "country",
356            RegionType::Island => "island",
357        };
358        println!(
359            "Region(48.8,2.3): {} ({})  area={:.0}km²",
360            r.name, rtype, r.area_km2
361        );
362    }
363
364    let elev_prov = ElevationProvider::global(30.0);
365    let elev_everest = elev_prov.sample(27.988, 86.925);
366    let notables = ElevationProvider::notable_elevations();
367    println!(
368        "Elevation(Everest)={:.0}m  notables={}",
369        elev_everest,
370        notables.len()
371    );
372
373    let bathy = BathymetryData::global(60.0);
374    let depth_mariana = bathy.sample(11.35, 142.2);
375    let is_ocean = bathy.is_ocean(0.0, -30.0);
376    let stats = BathymetryData::ocean_stats();
377    let basins = OceanBasin::major_basins();
378    println!(
379        "Bathy(Mariana)={:.0}m  is_ocean(0,-30)={}  basins={}  avg_depth={:.0}m",
380        depth_mariana,
381        is_ocean,
382        basins.len(),
383        stats.avg_depth_m
384    );
385
386    let mut hmap = Heightmap::new(64);
387    hmap.generate_procedural(42, 6, 0.5, 2.0);
388    let h_sample = hmap.sample(45.0, 10.0);
389    let r_at = hmap.radius_at(45.0, 10.0);
390    println!(
391        "Heightmap(64): sample(45,10)={:.1}m  radius={:.0}m",
392        h_sample, r_at
393    );
394
395    let config = LodConfig::default();
396    let mut lod = LodTerrain::new(config);
397    lod.update([6_371_000.0 + 1000.0, 0.0, 0.0]);
398
399    let mesh = TerrainMesh::from_region(0.0, 10.0, 0.0, 10.0, 8, &|lat, lon| hmap.sample(lat, lon));
400    println!(
401        "Mesh: {} vertices  {} triangles",
402        mesh.vertex_count(),
403        mesh.triangle_count()
404    );
405
406    let classifier = BiomeClassifier::default();
407    let biome = classifier.classify(2000.0, 45.0, 0.6);
408    let biome_name = match biome {
409        Biome::Ocean => "ocean",
410        Biome::Beach => "beach",
411        Biome::Desert => "desert",
412        Biome::Grassland => "grassland",
413        Biome::Forest => "forest",
414        Biome::Tundra => "tundra",
415        Biome::Snow => "snow",
416        Biome::Mountain => "mountain",
417        Biome::Volcanic => "volcanic",
418        Biome::Taiga => "taiga",
419    };
420    let splat = classifier.splat(2000.0, 45.0, 0.6);
421    println!(
422        "Biome(2000m,45°,0.6)={}  splat_top={:.2}",
423        biome_name, splat.weights[0].1
424    );
425
426    let face = Face::PosZ;
427    println!("LOD face={:?}", face);
428
429    let mat_ocean = PbrMaterial::ocean();
430    let mat_grass = PbrMaterial::grassland();
431    let mat_desert = PbrMaterial::desert();
432    let mat_snow = PbrMaterial::snow();
433    let mat_rock = PbrMaterial::rock();
434    let mat_volc = PbrMaterial::volcanic();
435    let mat_forest = PbrMaterial::forest();
436    let mat_ice = PbrMaterial::ice();
437    println!(
438        "Materials: ocean_rough={:.2}  grass_rough={:.2}  desert_metal={:.2}  \
439        snow_rough={:.2}  rock={:.2}  volc={:.2}  forest={:.2}  ice={:.2}",
440        mat_ocean.roughness,
441        mat_grass.roughness,
442        mat_desert.metallic,
443        mat_snow.roughness,
444        mat_rock.roughness,
445        mat_volc.roughness,
446        mat_forest.roughness,
447        mat_ice.roughness
448    );
449
450    let sh_terrain = ShaderData::terrain();
451    let sh_atm = ShaderData::atmosphere();
452    let sh_ocean = ShaderData::ocean();
453    println!(
454        "Shaders: terrain={} uniforms  atm={}  ocean={}",
455        sh_terrain.uniforms.len(),
456        sh_atm.uniforms.len(),
457        sh_ocean.uniforms.len()
458    );
459    for u in &sh_terrain.uniforms {
460        let val = match &u.value {
461            UniformValue::Float(f) => format!("f{:.2}", f),
462            UniformValue::Vec3(v) => format!("v3({:.1},{:.1},{:.1})", v[0], v[1], v[2]),
463            UniformValue::Vec4(v) => format!("v4({:.1},{:.1},{:.1},{:.1})", v[0], v[1], v[2], v[3]),
464            UniformValue::Mat4(_) => "mat4".to_string(),
465            UniformValue::Int(i) => format!("i{}", i),
466        };
467        println!("  {}={}", u.name, val);
468    }
469
470    let atm_params = AtmosphereParams::default();
471    let ray_dens = atm_params.rayleigh_density(4000.0);
472    let mie_dens = atm_params.mie_density(600.0);
473    let scat_r = atm_params.scatter_rayleigh(0.5);
474    let scat_m = atm_params.scatter_mie(0.5);
475    let sky = atm_params.sky_color([0.0, 1.0, 0.0], [0.0, 0.5, 0.866], 0.0);
476    println!(
477        "Atmosphere: ρ_ray(4km)={:.4}  ρ_mie(600m)={:.4}  scat_r={:.4}  scat_m={:.4}",
478        ray_dens, mie_dens, scat_r, scat_m
479    );
480    println!("  sky=({:.4},{:.4},{:.4})", sky[0], sky[1], sky[2]);
481
482    let ocean_rdr = OceanParams::default();
483    let phillips = ocean_rdr.phillips_spectrum(0.1, 0.05);
484    let spectrum = ocean_rdr.generate_spectrum();
485    let disp = ocean_rdr.dispersion(0.1);
486    println!(
487        "Ocean render: phillips(0.1,0.05)={:.6}  spectrum_size={}  disp={:.4}",
488        phillips,
489        spectrum.heights.len(),
490        disp
491    );
492
493    let cumulus = CloudLayer::cumulus();
494    let stratus = CloudLayer::stratus();
495    let cirrus = CloudLayer::cirrus();
496    let cb = CloudLayer::cumulonimbus();
497    let cloud_type = |ct: &CloudType| match ct {
498        CloudType::Cumulus => "Cu",
499        CloudType::Stratus => "St",
500        CloudType::Cirrus => "Ci",
501        CloudType::Cumulonimbus => "Cb",
502        CloudType::Altocumulus => "Ac",
503    };
504    println!(
505        "Clouds: {} base={:.0}m  {} base={:.0}m  {} base={:.0}m  {} base={:.0}m",
506        cloud_type(&cumulus.cloud_type),
507        cumulus.base_altitude_m,
508        cloud_type(&stratus.cloud_type),
509        stratus.base_altitude_m,
510        cloud_type(&cirrus.cloud_type),
511        cirrus.base_altitude_m,
512        cloud_type(&cb.cloud_type),
513        cb.base_altitude_m
514    );
515
516    let mut clouds = CloudSystem::earth_default();
517    clouds.step(3600.0);
518    let density = clouds.sample_density(5000.0, 45.0, 10.0);
519    println!("CloudSystem: density@5km={:.4}", density);
520
521    let rainforest = tropical_rainforest();
522    let boreal = boreal_forest();
523    let reef = coral_reef();
524    println!("Ecosystems:");
525    for eco in [&rainforest, &boreal, &reef] {
526        let sp = eco.total_species();
527        let ind = eco.total_individuals();
528        let shannon = eco.shannon_index();
529        let simpson = eco.simpson_diversity();
530        let area_sp = eco.expected_species_from_area(0.25, 100.0);
531        let npp = eco.total_npp_gc_yr();
532        let turnover = eco.biomass_turnover_time_yr(15_000.0);
533        println!(
534            "  {} species  {} ind  H'={:.2}  D={:.4}  SAR={:.0}  NPP={:.0}  τ={:.1}yr",
535            sp, ind, shannon, simpson, area_sp, npp, turnover
536        );
537    }
538
539    let tb = tropical_broadleaf();
540    let tg = temperate_grassland();
541    let cf = coniferous_forest();
542    let photo = tb.photosynthesis_rate(400.0, 1500.0, 28.0);
543    let canopy = tb.canopy_photosynthesis(photo);
544    let transp = tg.transpiration_mm_day(1.5, 0.3);
545    let npp_veg = cf.npp_kgc_m2_yr(15.0);
546    let crt = cf.carbon_residence_time_yr(npp_veg);
547    let light = beer_lambert_light_extinction(2000.0, 4.0, 0.5);
548    println!(
549        "Vegetation: photo={:.2}  canopy={:.2}  transp={:.2}mm/d  NPP={:.3}  CRT={:.1}yr  light={:.1}",
550        photo, canopy, transp, npp_veg, crt, light
551    );
552
553    let elephant = african_elephant();
554    let whale = blue_whale();
555    let gr_el = elephant.growth_rate();
556    let proj = elephant.project_forward(10.0);
557    let metab = whale.metabolic_rate_w();
558    let range = elephant.home_range_km2();
559    let gen_time = whale.generation_time_years();
560    let lifespan = elephant.max_lifespan_years();
561    println!(
562        "Elephant: r={:.3}  N(+10yr)={:.0}  range={:.1}km²  lifespan={:.0}yr",
563        gr_el, proj, range, lifespan
564    );
565    println!("Whale: metab={:.0}W  gen={:.1}yr", metab, gen_time);
566
567    let mut pp = PredatorPrey {
568        prey: african_elephant(),
569        predator: blue_whale(),
570        attack_rate: 0.01,
571        conversion_efficiency: 0.1,
572        predator_death_rate: 0.05,
573    };
574    let prey_r = pp.prey_growth_rate();
575    let pred_r = pp.predator_growth_rate();
576    pp.step(0.1);
577    println!(
578        "Lotka-Volterra: prey_r={:.1}  pred_r={:.1}  -> prey={:.1}  pred={:.2}",
579        prey_r, pred_r, pp.prey.count, pp.predator.count
580    );
581}

Source

pub fn sample_density(&self, altitude_m: f64, _lat: f64, _lon: f64) -> f64

Examples found in repository ?

examples/tidal_sim.rs (line 518)

47fn main() {
48    let orbit = EarthOrbit::new();
49    let period_days = orbit.orbital_period_days();
50    let v_perihelion = orbit.velocity_at_distance(orbit.perihelion_m());
51    let v_aphelion = orbit.velocity_at_distance(orbit.aphelion_m());
52    let energy = orbit.specific_orbital_energy();
53    let ang_mom = orbit.specific_angular_momentum();
54    let v_esc = EarthOrbit::escape_velocity_at_surface();
55    let f_sun = orbit.gravitational_force_sun();
56    let r = orbit.current_radius();
57    let v_mean = orbit.mean_orbital_velocity();
58    println!(
59        "Orbit: P={:.2}d  v_peri={:.0}m/s  v_aph={:.0}m/s  v_mean={:.0}m/s",
60        period_days, v_perihelion, v_aphelion, v_mean
61    );
62    println!(
63        "  E={:.2e}J/kg  L={:.2e}m²/s  v_esc={:.0}m/s  F_sun={:.2e}N  r={:.3e}m",
64        energy, ang_mom, v_esc, f_sun, r
65    );
66    println!(
67        "  a={:.3e}m  e={}  i={}°  Ω={}°  ω={}°",
68        SEMI_MAJOR_AXIS,
69        ECCENTRICITY,
70        INCLINATION_DEG,
71        LONGITUDE_ASCENDING_NODE_DEG,
72        ARGUMENT_PERIHELION_DEG
73    );
74
75    let rot = EarthRotation::new();
76    let v_equator = rot.surface_velocity_at_latitude(0.0);
77    let v_45 = rot.surface_velocity_at_latitude(45.0);
78    let accel = rot.centripetal_acceleration_at_latitude(0.0);
79    let coriolis = rot.coriolis_parameter(45.0);
80    let moi = rot.moment_of_inertia();
81    let ke = rot.rotational_kinetic_energy();
82    let l_rot = rot.angular_momentum();
83    let prec = rot.precession_rate_rad_per_year();
84    let daylen_summer = rot.day_length_variation_due_to_tilt(172, 60.0);
85    let daylen_winter = rot.day_length_variation_due_to_tilt(355, 60.0);
86    println!(
87        "Rotation: v_eq={:.0}m/s  v_45={:.0}m/s  a_c={:.4}m/s²  f_45={:.5e}",
88        v_equator, v_45, accel, coriolis
89    );
90    println!(
91        "  I={:.3e}kg·m²  KE={:.3e}J  L={:.3e}kg·m²/s  prec={:.2e}rad/yr",
92        moi, ke, l_rot, prec
93    );
94    println!(
95        "  day@60°N: summer={:.1}h  winter={:.1}h",
96        daylen_summer, daylen_winter
97    );
98    println!(
99        "  sidereal={:.4}s  solar={}s  tilt={}°={:.4}rad  prec_period={}yr",
100        SIDEREAL_DAY_S, SOLAR_DAY_S, AXIAL_TILT_DEG, AXIAL_TILT_RAD, PRECESSION_PERIOD_YEARS
101    );
102
103    let moon_tide = TidalForce::from_moon();
104    let sun_tide = TidalForce::from_sun();
105    let a_moon = moon_tide.tidal_acceleration();
106    let a_sun = sun_tide.tidal_acceleration();
107    let pot = moon_tide.tidal_potential(0.0);
108    let bulge_moon = moon_tide.tidal_bulge_height();
109    let bulge_sun = sun_tide.tidal_bulge_height();
110    let grav = moon_tide.gravitational_attraction();
111    let spring = spring_tide_amplitude();
112    let neap = neap_tide_amplitude();
113    let ratio = lunar_to_solar_tide_ratio();
114    println!(
115        "Tides: a_moon={:.2e}  a_sun={:.2e}  ratio={:.2}",
116        a_moon, a_sun, ratio
117    );
118    println!(
119        "  bulge_moon={:.3}m  bulge_sun={:.3}m  spring={:.3}m  neap={:.3}m",
120        bulge_moon, bulge_sun, spring, neap
121    );
122    println!("  potential(0)={:.2e}  F_grav={:.2e}N", pot, grav);
123
124    let chix = chicxulub_equivalent();
125    let tung = tunguska_equivalent();
126    let ke_chix = chix.kinetic_energy_mt();
127    let crater = chix.crater_diameter_m(2700.0);
128    let fireball = chix.fireball_radius_m();
129    let ejecta = chix.ejecta_volume_m3(2700.0);
130    let v_imp = chix.impact_velocity();
131    let ke_tung = tung.kinetic_energy_mt();
132    println!(
133        "Chicxulub: {:.2e}Mt  crater={:.0}m  fireball={:.0}m  ejecta={:.2e}m³  v_imp={:.0}m/s",
134        ke_chix, crater, fireball, ejecta, v_imp
135    );
136    println!("Tunguska: {:.2e}Mt", ke_tung);
137
138    let mut moon = MoonState::new();
139    let source = match moon.source.get() {
140        MoonSource::Binary => "binary",
141        MoonSource::Simulation => "simulated",
142    };
143    let pos0 = moon.position();
144    moon.step(3600.0);
145    let pos1 = moon.position();
146    let g_moon = moon.gravity_at(EARTH_MOON_DISTANCE);
147    println!(
148        "Moon({}): pos0=({:.0},{:.0})  pos1h=({:.0},{:.0})  g@dist={:.4e}m/s²",
149        source, pos0.0, pos0.1, pos1.0, pos1.1, g_moon
150    );
151    println!(
152        "  M={:.3e}kg  R={:.4e}m  d={:.3e}m  P={:.0}s",
153        LUNAR_MASS, LUNAR_RADIUS, EARTH_MOON_DISTANCE, LUNAR_ORBITAL_PERIOD
154    );
155
156    let iss = ArtificialSatellite::leo("ISS", 420_000.0, 408_000.0);
157    let geo = ArtificialSatellite::geo("GEO-Sat", 300.0);
158    let custom = ArtificialSatellite::new("Molniya", 1500.0, 500_000.0, 0.7, 1.1);
159    let p_iss = iss.orbital_period_s();
160    let v_iss = iss.orbital_velocity_ms();
161    let p_geo = geo.orbital_period_s();
162    let g_surf = iss.gravity_at_surface();
163    let pos_iss = iss.position();
164    let pos_geo = geo.position();
165
166    let orbit_type = |sat: &ArtificialSatellite| -> &str {
167        match sat.orbit_type {
168            OrbitType::LEO => "LEO",
169            OrbitType::MEO => "MEO",
170            OrbitType::GEO => "GEO",
171            OrbitType::HEO => "HEO",
172            OrbitType::Custom => "Custom",
173        }
174    };
175
176    println!(
177        "ISS({}): P={:.0}s  v={:.0}m/s  g_surf={:.2}  pos=({:.0},{:.0},{:.0})",
178        orbit_type(&iss),
179        p_iss,
180        v_iss,
181        g_surf,
182        pos_iss.0,
183        pos_iss.1,
184        pos_iss.2
185    );
186    println!(
187        "GEO({}): P={:.0}s  pos=({:.0},{:.0},{:.0})",
188        orbit_type(&geo),
189        p_geo,
190        pos_geo.0,
191        pos_geo.1,
192        pos_geo.2
193    );
194    println!(
195        "Custom({}): e={:.1}  i={:.1}rad",
196        orbit_type(&custom),
197        custom.eccentricity,
198        custom.inclination_rad
199    );
200
201    let mut constellation = Constellation::new("Starlink-test");
202    constellation.add(ArtificialSatellite::leo("S1", 260.0, 550_000.0));
203    constellation.add(ArtificialSatellite::leo("S2", 260.0, 550_000.0));
204    let mut sat_stepped = iss;
205    sat_stepped.step(60.0);
206    constellation.step_all(60.0);
207    let positions = constellation.positions();
208    println!(
209        "Constellation: {} sats  positions={}",
210        positions.len(),
211        positions.len()
212    );
213
214    let dt = DateTime::new(2024, 6, 21, 12, 0, 0.0);
215    let jd = dt.to_julian_date();
216    let dt_back = DateTime::from_julian_date(jd);
217    let unix = dt.to_unix_timestamp();
218    let dt_unix = DateTime::from_unix_timestamp(unix);
219    println!(
220        "DateTime: {}-{}-{} {}:{}:{:.0}  JD={:.4}  unix={:.0}",
221        dt.year, dt.month, dt.day, dt.hour, dt.minute, dt.second, jd, unix
222    );
223    println!(
224        "  roundtrip JD: {}-{}-{}  roundtrip unix: {}-{}-{}",
225        dt_back.year, dt_back.month, dt_back.day, dt_unix.year, dt_unix.month, dt_unix.day
226    );
227    println!(
228        "  J2000_JD={} UNIX_JD={} SPD={}",
229        J2000_JD, UNIX_EPOCH_JD, SECONDS_PER_DAY
230    );
231
232    let mut epoch = Epoch::j2000();
233    let epoch_mjd = Epoch::from_mjd(51544.5);
234    let centuries = epoch.centuries_since_j2000();
235    let days = epoch.days_since_j2000();
236    let gmst = epoch.gmst_degrees();
237    let mjd = epoch.to_mjd();
238    epoch.advance_days(365.25);
239    let centuries_1yr = epoch.centuries_since_j2000();
240    let mut epoch2 = Epoch::from_jd(J2000_JD);
241    epoch2.advance_seconds(86400.0);
242    println!(
243        "Epoch: T0={:.4}  days={:.1}  GMST={:.2}°  MJD={:.1}",
244        centuries, days, gmst, mjd
245    );
246    println!(
247        "  +1yr: T={:.6}  from_mjd.jd={:.1}  epoch2.days={:.1}",
248        centuries_1yr,
249        epoch_mjd.julian_date,
250        epoch2.days_since_j2000()
251    );
252    println!(
253        "  J2000={} J1950={} MJD_OFF={}",
254        J2000_EPOCH, J1950_EPOCH, MJD_OFFSET
255    );
256
257    let mut ts = TimeScale::realtime();
258    ts.step(1.0);
259    let sim_dt = ts.simulation_dt(1.0);
260    let hours = ts.sim_hours();
261    let days_ts = ts.sim_days();
262    let years = ts.sim_years();
263    ts.pause();
264    ts.resume();
265    ts.toggle_pause();
266    ts.set_speed(100.0);
267    let mut ff = TimeScale::fast_forward(10.0);
268    ff.step(1.0);
269    let sm = TimeScale::slow_motion(2.0);
270    println!(
271        "TimeScale: dt={:.1}  h={:.6}  d={:.8}  yr={:.10}  ff_sim={:.1}s  sm_speed={:.1}",
272        sim_dt, hours, days_ts, years, ff.simulation_time_s, sm.speed_multiplier
273    );
274
275    let sun_pos = SolarPosition::compute(jd, 48.8, 2.3);
276    let above = sun_pos.is_above_horizon();
277    let dist_sun = sun_pos.distance_m();
278    println!(
279        "Solar: el={:.1}°  az={:.1}°  above={}  dist={:.3e}m  tilt={}°",
280        sun_pos.elevation_deg, sun_pos.azimuth_deg, above, dist_sun, EARTH_AXIAL_TILT_DEG
281    );
282
283    let dnc = DayNightCycle::new(jd);
284    let state_paris = dnc.state_at(48.8, 2.3);
285    let state_text = match state_paris {
286        DaylightState::Day => "day",
287        DaylightState::Night => "night",
288        DaylightState::CivilTwilight => "civil_twilight",
289        DaylightState::NauticalTwilight => "nautical_twilight",
290        DaylightState::AstronomicalTwilight => "astro_twilight",
291    };
292    let terminator = dnc.terminator_points(36);
293    let ambient = dnc.ambient_light(48.8, 2.3);
294    println!(
295        "DayNight: Paris={}  ambient={:.2}  terminator_pts={}",
296        state_text,
297        ambient,
298        terminator.len()
299    );
300
301    let season_state = season_at(jd, 48.8);
302    let season_name = match season_state.season_north {
303        Season::Spring => "spring",
304        Season::Summer => "summer",
305        Season::Autumn => "autumn",
306        Season::Winter => "winter",
307    };
308    let south_name = match season_state.season_south {
309        Season::Spring => "spring",
310        Season::Summer => "summer",
311        Season::Autumn => "autumn",
312        Season::Winter => "winter",
313    };
314    let (sub_lat, sub_lon) = subsolar_point(jd);
315    println!(
316        "Season: N={}  S={}  decl={:.2}°  day_h={:.2}  subsolar=({:.1},{:.1})",
317        season_name,
318        south_name,
319        season_state.solar_declination_deg,
320        season_state.day_length_hours,
321        sub_lat,
322        sub_lon
323    );
324    println!(
325        "  tilt={}°  trop_year={}d  vernal_jd={}",
326        SEASONS_TILT, TROPICAL_YEAR_DAYS, VERNAL_EQUINOX_JD
327    );
328
329    let paris = LatLon::new(48.8566, 2.3522, 35.0);
330    let tokyo = LatLon::new(35.6762, 139.6503, 40.0);
331    let ecef_paris = paris.to_ecef();
332    let cart = paris.to_cartesian_simple();
333    let dist_pt = paris.distance_to(&tokyo);
334    let latlon_back = ecef_paris.to_latlon();
335    println!(
336        "Paris->ECEF: ({:.0},{:.0},{:.0})  cart=({:.0},{:.0},{:.0})",
337        ecef_paris.x, ecef_paris.y, ecef_paris.z, cart[0], cart[1], cart[2]
338    );
339    println!(
340        "  dist(Paris-Tokyo)={:.0}m  roundtrip_lat={:.4}°",
341        dist_pt, latlon_back.lat_deg
342    );
343    println!(
344        "  f={:.9}  a={:.0}m  b={:.3}m",
345        EARTH_FLATTENING, EARTH_SEMI_MAJOR_M, EARTH_SEMI_MINOR_M
346    );
347
348    let regions = RegionDatabase::continents();
349    let region = regions.point_in_region(48.8, 2.3);
350    if let Some(r) = region {
351        let rtype = match r.region_type {
352            RegionType::Continent => "continent",
353            RegionType::Ocean => "ocean",
354            RegionType::Sea => "sea",
355            RegionType::Country => "country",
356            RegionType::Island => "island",
357        };
358        println!(
359            "Region(48.8,2.3): {} ({})  area={:.0}km²",
360            r.name, rtype, r.area_km2
361        );
362    }
363
364    let elev_prov = ElevationProvider::global(30.0);
365    let elev_everest = elev_prov.sample(27.988, 86.925);
366    let notables = ElevationProvider::notable_elevations();
367    println!(
368        "Elevation(Everest)={:.0}m  notables={}",
369        elev_everest,
370        notables.len()
371    );
372
373    let bathy = BathymetryData::global(60.0);
374    let depth_mariana = bathy.sample(11.35, 142.2);
375    let is_ocean = bathy.is_ocean(0.0, -30.0);
376    let stats = BathymetryData::ocean_stats();
377    let basins = OceanBasin::major_basins();
378    println!(
379        "Bathy(Mariana)={:.0}m  is_ocean(0,-30)={}  basins={}  avg_depth={:.0}m",
380        depth_mariana,
381        is_ocean,
382        basins.len(),
383        stats.avg_depth_m
384    );
385
386    let mut hmap = Heightmap::new(64);
387    hmap.generate_procedural(42, 6, 0.5, 2.0);
388    let h_sample = hmap.sample(45.0, 10.0);
389    let r_at = hmap.radius_at(45.0, 10.0);
390    println!(
391        "Heightmap(64): sample(45,10)={:.1}m  radius={:.0}m",
392        h_sample, r_at
393    );
394
395    let config = LodConfig::default();
396    let mut lod = LodTerrain::new(config);
397    lod.update([6_371_000.0 + 1000.0, 0.0, 0.0]);
398
399    let mesh = TerrainMesh::from_region(0.0, 10.0, 0.0, 10.0, 8, &|lat, lon| hmap.sample(lat, lon));
400    println!(
401        "Mesh: {} vertices  {} triangles",
402        mesh.vertex_count(),
403        mesh.triangle_count()
404    );
405
406    let classifier = BiomeClassifier::default();
407    let biome = classifier.classify(2000.0, 45.0, 0.6);
408    let biome_name = match biome {
409        Biome::Ocean => "ocean",
410        Biome::Beach => "beach",
411        Biome::Desert => "desert",
412        Biome::Grassland => "grassland",
413        Biome::Forest => "forest",
414        Biome::Tundra => "tundra",
415        Biome::Snow => "snow",
416        Biome::Mountain => "mountain",
417        Biome::Volcanic => "volcanic",
418        Biome::Taiga => "taiga",
419    };
420    let splat = classifier.splat(2000.0, 45.0, 0.6);
421    println!(
422        "Biome(2000m,45°,0.6)={}  splat_top={:.2}",
423        biome_name, splat.weights[0].1
424    );
425
426    let face = Face::PosZ;
427    println!("LOD face={:?}", face);
428
429    let mat_ocean = PbrMaterial::ocean();
430    let mat_grass = PbrMaterial::grassland();
431    let mat_desert = PbrMaterial::desert();
432    let mat_snow = PbrMaterial::snow();
433    let mat_rock = PbrMaterial::rock();
434    let mat_volc = PbrMaterial::volcanic();
435    let mat_forest = PbrMaterial::forest();
436    let mat_ice = PbrMaterial::ice();
437    println!(
438        "Materials: ocean_rough={:.2}  grass_rough={:.2}  desert_metal={:.2}  \
439        snow_rough={:.2}  rock={:.2}  volc={:.2}  forest={:.2}  ice={:.2}",
440        mat_ocean.roughness,
441        mat_grass.roughness,
442        mat_desert.metallic,
443        mat_snow.roughness,
444        mat_rock.roughness,
445        mat_volc.roughness,
446        mat_forest.roughness,
447        mat_ice.roughness
448    );
449
450    let sh_terrain = ShaderData::terrain();
451    let sh_atm = ShaderData::atmosphere();
452    let sh_ocean = ShaderData::ocean();
453    println!(
454        "Shaders: terrain={} uniforms  atm={}  ocean={}",
455        sh_terrain.uniforms.len(),
456        sh_atm.uniforms.len(),
457        sh_ocean.uniforms.len()
458    );
459    for u in &sh_terrain.uniforms {
460        let val = match &u.value {
461            UniformValue::Float(f) => format!("f{:.2}", f),
462            UniformValue::Vec3(v) => format!("v3({:.1},{:.1},{:.1})", v[0], v[1], v[2]),
463            UniformValue::Vec4(v) => format!("v4({:.1},{:.1},{:.1},{:.1})", v[0], v[1], v[2], v[3]),
464            UniformValue::Mat4(_) => "mat4".to_string(),
465            UniformValue::Int(i) => format!("i{}", i),
466        };
467        println!("  {}={}", u.name, val);
468    }
469
470    let atm_params = AtmosphereParams::default();
471    let ray_dens = atm_params.rayleigh_density(4000.0);
472    let mie_dens = atm_params.mie_density(600.0);
473    let scat_r = atm_params.scatter_rayleigh(0.5);
474    let scat_m = atm_params.scatter_mie(0.5);
475    let sky = atm_params.sky_color([0.0, 1.0, 0.0], [0.0, 0.5, 0.866], 0.0);
476    println!(
477        "Atmosphere: ρ_ray(4km)={:.4}  ρ_mie(600m)={:.4}  scat_r={:.4}  scat_m={:.4}",
478        ray_dens, mie_dens, scat_r, scat_m
479    );
480    println!("  sky=({:.4},{:.4},{:.4})", sky[0], sky[1], sky[2]);
481
482    let ocean_rdr = OceanParams::default();
483    let phillips = ocean_rdr.phillips_spectrum(0.1, 0.05);
484    let spectrum = ocean_rdr.generate_spectrum();
485    let disp = ocean_rdr.dispersion(0.1);
486    println!(
487        "Ocean render: phillips(0.1,0.05)={:.6}  spectrum_size={}  disp={:.4}",
488        phillips,
489        spectrum.heights.len(),
490        disp
491    );
492
493    let cumulus = CloudLayer::cumulus();
494    let stratus = CloudLayer::stratus();
495    let cirrus = CloudLayer::cirrus();
496    let cb = CloudLayer::cumulonimbus();
497    let cloud_type = |ct: &CloudType| match ct {
498        CloudType::Cumulus => "Cu",
499        CloudType::Stratus => "St",
500        CloudType::Cirrus => "Ci",
501        CloudType::Cumulonimbus => "Cb",
502        CloudType::Altocumulus => "Ac",
503    };
504    println!(
505        "Clouds: {} base={:.0}m  {} base={:.0}m  {} base={:.0}m  {} base={:.0}m",
506        cloud_type(&cumulus.cloud_type),
507        cumulus.base_altitude_m,
508        cloud_type(&stratus.cloud_type),
509        stratus.base_altitude_m,
510        cloud_type(&cirrus.cloud_type),
511        cirrus.base_altitude_m,
512        cloud_type(&cb.cloud_type),
513        cb.base_altitude_m
514    );
515
516    let mut clouds = CloudSystem::earth_default();
517    clouds.step(3600.0);
518    let density = clouds.sample_density(5000.0, 45.0, 10.0);
519    println!("CloudSystem: density@5km={:.4}", density);
520
521    let rainforest = tropical_rainforest();
522    let boreal = boreal_forest();
523    let reef = coral_reef();
524    println!("Ecosystems:");
525    for eco in [&rainforest, &boreal, &reef] {
526        let sp = eco.total_species();
527        let ind = eco.total_individuals();
528        let shannon = eco.shannon_index();
529        let simpson = eco.simpson_diversity();
530        let area_sp = eco.expected_species_from_area(0.25, 100.0);
531        let npp = eco.total_npp_gc_yr();
532        let turnover = eco.biomass_turnover_time_yr(15_000.0);
533        println!(
534            "  {} species  {} ind  H'={:.2}  D={:.4}  SAR={:.0}  NPP={:.0}  τ={:.1}yr",
535            sp, ind, shannon, simpson, area_sp, npp, turnover
536        );
537    }
538
539    let tb = tropical_broadleaf();
540    let tg = temperate_grassland();
541    let cf = coniferous_forest();
542    let photo = tb.photosynthesis_rate(400.0, 1500.0, 28.0);
543    let canopy = tb.canopy_photosynthesis(photo);
544    let transp = tg.transpiration_mm_day(1.5, 0.3);
545    let npp_veg = cf.npp_kgc_m2_yr(15.0);
546    let crt = cf.carbon_residence_time_yr(npp_veg);
547    let light = beer_lambert_light_extinction(2000.0, 4.0, 0.5);
548    println!(
549        "Vegetation: photo={:.2}  canopy={:.2}  transp={:.2}mm/d  NPP={:.3}  CRT={:.1}yr  light={:.1}",
550        photo, canopy, transp, npp_veg, crt, light
551    );
552
553    let elephant = african_elephant();
554    let whale = blue_whale();
555    let gr_el = elephant.growth_rate();
556    let proj = elephant.project_forward(10.0);
557    let metab = whale.metabolic_rate_w();
558    let range = elephant.home_range_km2();
559    let gen_time = whale.generation_time_years();
560    let lifespan = elephant.max_lifespan_years();
561    println!(
562        "Elephant: r={:.3}  N(+10yr)={:.0}  range={:.1}km²  lifespan={:.0}yr",
563        gr_el, proj, range, lifespan
564    );
565    println!("Whale: metab={:.0}W  gen={:.1}yr", metab, gen_time);
566
567    let mut pp = PredatorPrey {
568        prey: african_elephant(),
569        predator: blue_whale(),
570        attack_rate: 0.01,
571        conversion_efficiency: 0.1,
572        predator_death_rate: 0.05,
573    };
574    let prey_r = pp.prey_growth_rate();
575    let pred_r = pp.predator_growth_rate();
576    pp.step(0.1);
577    println!(
578        "Lotka-Volterra: prey_r={:.1}  pred_r={:.1}  -> prey={:.1}  pred={:.2}",
579        prey_r, pred_r, pp.prey.count, pp.predator.count
580    );
581}