We simulated hadronic showers at EeV energies with CORSIKA to study the massive quark production and decay modes. In this regard, we specifically studied charm meson production and their decay channels with SYBILL model. This understanding would possibly help implement bottom quark production and their decay in these models. The massive quark decays contribute to the atmospheric high-energy gamma-rays and neutrinos. Hence, this study would be essential to significantly observe high-energy gamma rays and neutrinos fluxes from astrophysical objects.